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We develop the formalism required to study the nonlinear interaction of modes in rotating Newto- 
nian stars, assuming that the mode amplitudes are only mildly nonlinear. The formalism is simpler 
than previous treatments of mode-mode interactions for spherical stars, and simplifies and corrects 
previous treatments for rotating stars. At linear order, we elucidate and extend slightly a formalism 
due to Schutz, show how to decompose a general motion of a rotating star into a sum over modes, 
and obtain uncoupled equations of motion for the mode amplitudes under the influence of an ex- 
ternal force. Nonlinear effects are added perturbatively via three-mode couplings, which suffices for 
moderate amplitude modal excitations; the formalism is easy to extend to higher order couplings. 
We describe a new, efficient way to compute the modal coupling coefficients, to zeroth order in the 
stellar rotation rate, using spin-weighted spherical harmonics. The formalism is general enough to 
allow computation of the initial trends in the evolution of the spin frequency and differential rotation 
of the background star. 

We apply this formalism to derive some properties of the coupling coefficients relevant to the 
nonlinear interactions of unstable r-modes in neutron stars, postponing numerical integrations of 
the coupled equations of motion to a later paper. First, we clarify some aspects of the expansion 
in stellar rotation frequency fl that is often used to compute approximate mode functions. We 
show that in zero-buoyancy stars, the rotational modes (those modes whose frequencies vanish as 
Q — > 0) are orthogonal to zeroth order in fi. From an astrophysical viewpoint, the most interesting 
result of this paper is that many couplings of r— modes to other rotational modes are small: either 
they vanish altogether because of various selection rules, or they vanish to lowest order in Q or in 
compressibility. In particular, in zero-buoyancy stars, the coupling of three r— modes is forbidden 
entirely and the coupling of two r-modes to one hybrid, or r-g rotational mode vanishes to zeroth 
order in rotation frequency. The coupling of any three rotational modes vanishes to zeroth order in 
compressibility and in £1. In nonzero-buoyancy stars, coupling of the r-modes to each other vanishes 
to zeroth order in Q, . Couplings to regular modes (those modes whose frequencies are finite in the 
limit il — > 0), such as /—modes, are not zero, but since the natural frequencies of these modes are 
relatively large in the slow rotation limit compared to those of the r-modes, energy transfer to those 
modes is not expected to be efficient. 



I. INTRODUCTION 



Various authors have proposed that r-modes might be linearly unstable inside fairly rapidly rotating neutron stars, 
depending on a competition between gravitational radiation reaction, which drives the instability, and viscous effects, 
which inhibit it; several other investigators have explored implications for the spin evolution of rapidly rotating new- 
born and accreting neutron stars fl|-|l4f|. An open question is how the instability saturates, and at what perturbation 
amplitude. This question is not only important from a theoretical standpoint, but also must be answered in order 
to assess whether or not the gravitational radiation emitted during the development of the instability is detectable. 
The saturation amplitude also determines the final spin rate of the neutron star when it exits the instability region of 
parameter space. While turbulence generated in the strong shear layer at the interface between the fluid core and solid 
crust has been shown to saturate the r-mode p"5[ , a crust will not be present initially for newly born neutron stars, 
and other hydrodynamical saturation mechanisms are still being investigated and may in fact be more important. 

There have been, to date, two fully nonlinear numerical calculations of the development of the instability, but 
neither completely settles the question of saturation, principally because there are practical limitations on what it 
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is feasible to compute. The first calculation, by Stergioulas and Font fuH , considered the evolution of the fluid 
inside a neutron star with no buoyancy force and fixed spacetime geometry. Their calculations began with a large 
amplitude excitation in the mode expected to be most unstable, the I = m — 2 r-mode, and followed the subsequent 
hydrodynamics for about 20 stellar rotation periods, approximately 25 ms for their calculations. During this time 
no substantial change in the amplitude of the modal excitation was seen, and no evidence for significant nonlinear 
excitation of other stellar normal modes was detected. Lindblom, Tohline and Vallisneri |]l7| subsequently performed 
simulations based on the Newtonian equations of hydrodynamics and gravitation, including an artificially enhanced 
gravitational radiation force consisting of the correct force jl8|,[l9| multiplied by a factor of about 4500 pofl . The 
linear instability growth timescale in their simulations was about f 3 rotation periods, rather than the correct value of 
> 10 4 rotation periods. The initial r— mode amplitude chosen by Lindblom, Tohline and Vallisneri was small (about 
O.f ), and, as expected, grew via the radiation-driven instability to large amplitude, where shocks ultimately formed, 
accompanied by a decrease in r— mode amplitude. Both calculations show that nonlinear mode-mode coupling is 
comparatively weak for r— modes, producing little effect even at large modal amplitude on timescales ~ f — 20 stellar 
rotation periods. 

A second, more analytical approach to this problem is possible. One assumes, at the outset, that the modes develop 
only moderately nonlinear amplitudes, and nonlinear effects are modeled perturbatively, via n— mode couplings. When 
amplitudes are small, relatively low order couplings, such as three-mode interactions, suffice, and (when dissipation 
and radiation reaction are ignored) the equations of motion for the fluid can be modeled as an infinite-dimensional 
Hamiltonian system with a polynomial interaction potential (which is cubic when only three-mode couplings are 
kept). This approach has the advantage of allowing explicit calculation of the strength of the lowest order coupling 
of a given r— mode to other modes of the star, but has the disadvantage of failing when mode amplitudes grow very 
large. Although it is possible that the r— mode instability saturates only in the fully nonlinear regime, neither of 
the numerical calculations establish this, because it is possible that modal couplings act on timescales longer than 
covered by the simulations, but shorter than the growth time (~ f 4 rotation periods) of the gravitational-radiation 
instability. The analytic approach we are taking should apply if the instability saturates before attaining very large 
amplitudes. Moreover, by computing the three-mode coupling coefficients for r— modes, we should be able to shed 
some light on why nonlinear coupling is so inefficient for them on dynamical timescales in the simulations. 

In this paper, we develop the tools needed to calculate the nonlinear evolution of the unstable modes in the weakly 
nonlinear regime. We concentrate on three-mode interactions for the most part, although some of the results we 
derive could also apply to interactions involving n > 3 modes of the star. 

We begin with a review of the linear perturbation theory in Sec. [n]. A principal goal of that section of the paper is 
to obtain equations of motion for modes of a rotating star acted on by an external force that are uncoupled from one 
another at linear order. The method we develop is based on a series of papers by Schutz and collaborators [pl|-p6[, 
but has not been applied consistently to r— modes before. The key results, for a reader wishing to skip the details, are 
the mode decomposition formulae (2.47) and ( |2.48| ) and the equation of motion for each mode (2.49). In Sec. Ill, we 
review the slow rotation approximation to mode functions and frequencies, with particular attention to the rotational 
modes, whose frequencies are zero in the limit of zero rotation frequency, f2. The key new result in this section is the 
fact that the rotational mode functions are orthogonal to leading order in the stellar angular velocity. 

In Sec. [V, we discuss the second order Lagrangian perturbation theory at the heart of our calculation of the evolution 
of mode s in t he weakly nonlinear regime. The key results of this section are (i) the coupled nonlinear equations of 
motion ( 4.17 ) for the mode coefficients, and (ii) the formula ( 4.2C ) for the three-mode coupling coefficients that can 
be evaluated once the modes have been found for a given background star model. In Sec. [v|, we discuss properties of 
the coupling coefficients. We show how general parity arguments and selection rules can be used to explain why some 
modes do not couple to one another, and also show that there are other couplings that are suppressed, in that they 
may vanish to lowest order in Q, or to lowest order in compressibility, but are nonzero more generally. In Sec. VI, we 



develop a new and efficient way of computing coupling coefficients that employs spin-weighted spherical harmonics 
p7t . The method greatly simplifies the calculation of the coupling coefficients for nonrotating stars as well as rotating 
stars, and should be useful for other astrophysical applications (see, e.g., Refs. [f28| , ^9|| ) . 

This paper concentrates on formal issues primarily, postponing numerical integration of the time-dependent evo- 
lution of a network of weakly interacting modes to a subsequent publication. From a formal viewpoint, the most 
substantial result of this paper is a formulation of the equations of motion that, at linear order, separates the re- 
sponses of individual modes to an external force, and, at the first nonlinear order, involves coupling coefficients that 
we can calculate relatively economically. From an astrophysical viewpoint, the most interesting result is that many 
interactions of r— modes with themselves can be shown to be either forbidden entirely by selection rules, or suppressed 
to lowest order in some parameter, such as rotation frequency or compressibility. By itself, this is not quite enough 
to explain the smallness of the coupling seen in numerical simulations, but we may already be seeing a hint of why 
the transfer of energy from r— modes to other modes is so inefficient on dynamical timescales. 
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II. LINEAR PERTURBATION THEORY FOR ROTATING STARS: THE SCHUTZ FORMALISM 



A. Overview 



In this section we develop the formalism we need in the linear Lagrangian perturbation theory of rotating stars. 
There are two aspects to the perturbation theory that will be essential for our purposes: 

(i) Solving for the frequencies uj a and mode functions £ Q (x)- This aspect is well understood and highly developed; 
see for example Ch. VI of Unno et al J|(|, as well as Refs. |pl|-^|. 

(ii) Decomposing general motions in a star into sums over modes, and obtaining equations of motion for the mode 
coefficients that are uncoupled at linear order when an external force acts on the star. 

The second aspect of the theory — deriving equations of motion for the mode coefficients — has not been well- 
understood in the past. The necessary theory is implicit in a series of papers by Schutz and collaborators [pll-p6|, 
but the formalism has not been applied widely in astrophysics. For example, a classic computation for which modal 
equations of motion would be useful is the excitation of the modes of a rotating star by the tidal field of a binary 
companion. While there have been many papers on this topic, we are not aware of any that have use d the correct 
uncoupled equations of motion. For example, Refs. |54j and |35| use a type of mode expansionjEq. ( 2.50 ) below] 
for which the equations of motion of the various mode amplitudes are coupled, and Refs. |3(J, |37J] and |3§|1 use the 
same mode expansion but simply drop the coupling terms between the different modes. Some authors have resorted 
to numerically solving the linearized hydrodynamic equations in the frequency domain, after factoring out a factor of 
e im<p—iut^ ra ther than using a mode expansion |39|-f4l|] . 

This section of the paper is devoted to presenting and explaining the equations of motion that result from Schutz 's 
formalism when one makes a number of simplifying assumptions that usually are valid in applications. We have 
attempted to make the presentation trans pare nt so as to be accessible to a wide audience. All the derivations are 
relegated to an appendix. We start in Sec. |IIB| by summarizing the governing equations. In Sec. II C we describe the 
form of the mode expansion and equations of motion for the mode coefficients. 



B. Governing equations 

We assume that the unperturbed background star is uniformly rotating with angular velocity SI. The fluid equations 
of motion in the corotating frame are 

^+V-(pu) = (2.1) 

and 

+ (u ■ V)u + 2SI x u + SI x (SI x x) = - + a ex t, (2.2) 

at p 

where p is the density, p the pressure, u the velocity and a oxt is any acceleration due to external forces. The 
gravitational potential <f> is given by V 2 = AirGp. In the background solution, u vanishes and p is time independent. 

Now consider linearized perturbations characterized by Eulerian perturbations Su, Sp and Sp. The first two of these 
are related to the linearized Lagrangian displacement £(x, t) by 

6p = -V> (p^) (2.3) 

and Su = £ (since u = 0). We assume that the perturbations are characterized by an adiabatic index T\, so that 

^ = r x ^, (2.4) 
p p 

where Ap = Sp + (£ • V)p and Ap = Sp + (£ • V)p are the Lagrangian perturbations of density and pressure^]. By 
combining these relations with linearized versions of the continuity and Euler equations (2.1) and (|2.2|), 



1 If the background star and the perturbations both obey the same one-parameter equation of state p = p(p), then Ti = 
(p/p)dp/dp depends only on the background density p. More generally, we can regard I\ = Ti(x) as a function of position 
determined by the structure of the backgroun d s tar. A very gen eral definition of the adiabatic index T 1 — I\(x) that is 
consistent with the phenomenological relation (2.4) is given in Eq. (J24) below. 
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5p + V • (p6u) = 0, (2.5) 

and 

Su + 20 x 5u = -Y^H + YPfip - S/Sd> + a oxt , (2.6) 
P P 

where V 2 5(f> = AirGSp, we obtain the equations of motion for linearized perturbations in the form [ fi2"| 

<j + B.£ + C-€ = a«t(x,t). (2.7) 

Here the operator B is given by 

B-£ = 2fJx£, (2.8) 

and the operator C is given by 

P (C • |) 4 = V,; I i/»V,;"*' :• + V lP V 3 e ~ Vj-p V 4 f 

+p^V j V i (j) + ^VjVi^rot + pV^</>, (2.9) 

where 

fc t (x) = -i(fixxf. (2.10) 

The C operator can be broken up into two pieces 

C = C a + C 6 (2.11) 
where C a is proportional to the vectorial Schwarzschild discriminant 

A-^-i?, (2.12, 



and C& is independent of A. We find 



and 



P Ti p 



C . C = vf^+^] (2.18-. 

V p p 



C,^ V fe.A) + *A + ^.A (2.14) 

V p J p p 

= ^SpTxA + -V (T lP $ ■ A) . (2.15) 
P P 

The buoyancy force described by the operator C& gives rise to g modes whose frequency scale is set by the Brunt- 
Vaisala frequency N (Brunt for short), defined by 

N 2 = ^-A. (2.16) 
P 

For much of this paper we will specialize to perturbations of stars without buoyancy, i.e., situations where A = 00. 
For example, if the background star and the perturbations obey one parameter equations of state with adiabatic 
indices T and Ti, respectively, then A oc T — I\, and so the zero-buoyancy case is T — T\. For ordinary stars, the 
Brunt frequency in the convection zone is approximately zero when the superadiabatic gradient is small. For neutron 



2 For such stars dp A dp = and so the background star satisfies p = p(p) for some equation of state. Furthermore the 
perturbations obey the same equation of state, so the star is barotropic. [A special case of this is a star of uniform specific 
entropy, i.e., an isentropic star.] For more general stars dp and dp need not be proportional. 
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stars, A = means a cold, zero-entropy gas in which the composition adjusts instantly to changes in density. In real 
neutron stars the adjustment speed is limited by beta and inverse-beta decays, which gives rise to g modes |43]j (see 
Appendix Q for a related discussion). Our analysis of r modes will make the zero- buoyancy approximation A = 0, 
which is expected to be good in the regime N <§C O. However, this inequality will only be marginally satisfied in 
newly born neutron stars, so it will be important for future analyses to analyze nonlinear couplings of r modes in 
nonzero-buoyancy stars. We cho ose to focus on the zero-buoyancy case for simplicity. 

The dynamical equation (2.7) forms the starting point for our discussions. We can find a large class of solutions 
for the case a cxt = of no forcing if we make the ansatz 



£(x,t) 



where ui is the rotating-frame frequency. This yields 

[-lu 2 - iuB + C] • £ = 0. 



(2.17) 



(2.18) 



The quadratic eigenvalue equation (2.18) is the standard equation that one solves to obtain the eigenfunctions £(x) 
and eigenfrequencies uj for rotating stars. A mode of the star will consist of a pair (£,u>) that satisfy Eq. ( |2.18|) . 



C. Mode decomposition formalism 

1. Nonrotating stars 

We start by recalling the standard mode decomposition formalism fo r non -rotating stars for which B = 0, in order 
to contrast it with the rotating case below. For non-rotating stars Eq. ( 2.18 ) reduces to 



(2.19) 



which represents a standard eigenvalue problem for to 2 . We define the inner product on the space Ji of complex vector 
functions £ = £(x) by 

(CO = / d 3 xp( X )^)* -r(x). (2.20) 
The operator C is Hermitian with respect to this inner product^]: for any elements £ and of H, 

<e,c-o = <c-e,o> ( 2 - 21 ) 

while B is anti-Hermitian 

(€,B-0 = -<B-€>0- (2-22) 
It follows that one can find a set {£ Q } °f eigenvectors of C with eigenvalues L0 2 a which are orthonormal 

(Sa > C/9> = <W (2.23) 
and which also form a basis for the vector space TL. Hence we can decompose any Lagrangian displacement £(x, t) as 



*(x,t)=X>a(*)€a(x), 



(2.24) 



3 Throughout this paper we assume that the space H. is finite dimensional, as it will be in any numerical computations, to 
avoid discussing subtleties relate d to infinite dimensional spaces. 



4 The sum over a in Eq. (2.24) is o ver a ll distinct solutions £ Q to Eq. (2. It), not over the larger set of distinct pairs 
as there is only one term in the sum ( 2.24 ) corresponding to the two solutions (£ a , u a ) and (£ a , — w a ). This will be important 
below. 
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where 



fe(*) = <€„>€(*)) 



(2.25) 



Finally, by inserting the expansion ( 2.24 ) into the dynamical equation (2/7) with B = and by contracting with £ a , 
one obtains the equation of motion 



q a {t) + uj 2 a q a {t) = (£ Q , a oxt (t)) 
which is just the forced harmonic oscillator equation. 



(2.26) 



2. Rotating stars 

Turn, now, to the corresponding formalism for rotating stars. We make a number of simplifying assumptions: 
• The background star is stable, so all frequencies are real. 



All the modes are non-degenerate: for each frequency ui the dimension of the space of eigenmodes is one. This 



assumption will be relaxed in Sec. II C 3 below 



A more general treatment, witho ut these assumptions, is given in Appendix 

The eigenvalue equation (2.18) can be reformulated as a standard eigenvalue equation [Eq. ( |A9[ ) below] for uj involv- 



ing a non-Hermitian operator. Non-Hermitian operators in general do not posess sufficiently many right eigenvectors 
to form a basis. However, for such operators there is a standard procedure described in linear algebra textbooks 
to obtain a basis by adding to the eigenvectors certain additional vectors, sometimes called generalized eigenvectors 
or Jordan chain vectors. Each actual eigenvector may have one or more generalized eigenvectors associated with it, 
forming a so-called Jordan chain. 

There are two classes of modes (£, ui) in a rotating star: modes associated with non-trivial Jordan chains (E5|, 
which we shall call Jordan-chain modes, and modes not associated with such Jordan chains (i.e., the length of the 
chain is zero; see Appendix ^ for more details). A mode that is just at the point of becoming unstable is always a 



Jordan-chain mode |25|. However, even for stable rotating stars, Jordan-chain modes of zero frequency are always 
present. One set of such modes corresponds to displacements of the star to nearby equilibria with slightly different 
angular velocities 

Q^n + 6Q(r±), (2.27) 

where r± is the distance from the rotation axis (see Appendix [j5] for a proof of this in the zero-buoyancy case)[]. 
Therefore, one cannot assume that no Jordan chains occur. In this section we shall nevertheless, for simplicity, 
describe the formalism that would apply if there were no Jordan-chain modes. This formalism should be useful, for 
example, in describing situations in which the Jordan-chain modes are unimportant dynamically. In the remainder of 
this paper we shall assume that the Jordan chain modes are unimportant for the process of saturation of r-modes^. 
A more complete treatment of the formalism, allowing Jordan-chain modes, is detailed in Appendix |a|. 

Given these assumptions, there are two obstacles to obtaining a mode d ecom position formalism similar to that for 



non-rotating stars. First, distinct mode s (£, u>) and satisfying Eq. ( 2.18 ) with uj ^ uj' need not be orthogonal 

with respect to the inner product ( |2.20|) . Although the operator 

L(w) = -uj 2 - iujB + C (2.28) 

is Hermitian for real uj |42[, the vectors £ and £' are eigenvectors of the two different operat ors L (uj) and L(w') and 



so need not be orthogonal. Nevertheless, it is possible to find a set {£, a } of solutions to Eq. ( 2.18 ) that form a basis 



of TC, so that every Lagrangian displacement can be uniquely decomposed as a superposition of the form 



5 There are also Jordan chain modes corresponding to tilting the axis of rotation of the star. 

6 However, note that several studies have found that the gravitational radiation reaction and/or hydrodynamic nonlinearities 
induce significant differential rotation in the star Ju^|l6l|l7t . A set of coupled equations for the growth of the r-mode as well as 
for the differential rotation can be formulated using the results of A ppen dix ^| Since the differential rotation corresponds to a 
set of Jordan chains of length one (Appendix^), Eqs. (A21) and (A.22) can be used to evolve the differential rotation in the 
linear regime; however the linear approximation breaks down after a time ~ 1/8^1. 
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€(x,t)=5Dg a (t)£ a (x) 



(2.29) 



The second obstacle is the following: There is in g enera l no choice of basis {£ Q } of eigenvectors for which the equations 
of motion of the coefficients q a (t) defined by Eq. (2.21) are uncoupled from one another. 

To circumvent this obstacle, it is necessary to use a phase space mode expansion instead of a configuration space 
mode e xpan sion, as pointed out by Dyson and Schutz pi) . Let's label the distinct non- Jordan-chain solutions^ 
to Eq. ( 2.18 ) as (£ a ,uja), so that 



The set of vectors 



- iuj A B + C] • £ A = 0. 



(2.30) 



(2.31) 



We deliberately use a different type of index here - 
as a reminder that the number of distinct A's is twice 



forms a basis for the space Tt © H of pairs of vectors [£, £'] 
capital Roman indices rather than lower case Greek indices - 
the dimension of 7i, whereas the number of distinct a's in the sums fl2.24j ) and Q2.29D is just the dimension of Ti.. 

The mode functions £ A are not orthogonal in general, that is, ($4 , £ B ) need not vanish when A B. We give an 
explicit example demonstrating the non-orthogonality in Sec. Q below. However, as shown by Ref. f22| , the modes 
(£^, w^) do obey a modified type of orthogonality relation, which is ^| 



for A^B. 

Now, at any time £, we can form the vector 



(wa+w b )<Ca.Cb) = 



€(*) 



(2.32) 



(2.33) 



from the Lagrangian displacement £(x, t) and its time derivative £(x, t). We can expand this vector in the basis (2.31) 



(2.34) 



Using the orthogonality relation (2.32), the expansion (2.34) can be inverted (see Appendix |A|) to obtain 



c A (t) = — {Za . + *€(*) + *B • c(t) 



where the constant 6a is given by 



(2.35) 



6,4 = (^,zB-^)+2^ 



(2.36) 



7 Here we define two solutions 
8 The left hand side of Eq. ([[.3 

£a( x ) ex P[ — WAt] and £ 2 ( x j£) = 
construction was given by Ref. 
9 Note that it follows from Eq 




uj) and to be distinct if lu 7^ a/ or if £ and £' are linearly independent, 

is proportional to the symplectic product W(£ 1 ,£ 2 ) defined by Ref. [ B 2 J , where ^(x, t) 
(x) exp[— iuat]. An alternative proof of the orthogonality relation ( |2,32| ) from a phase space 
which we review in Appendix | 
34h that 



(ca + IUaCa)£j 



0. 



However, it does not follow that ca + iloaca = for each A since the set {£ A } is not a basis for 7i; it has twice too many basis 
elements. 

10 Note that adjusting the normalization of £ A will change the magnitude but not the sign of b A . Hence we cannot find a basis 
that achieves 6 a = 1 for all A, since some of the 6a 's are negative and some are positive. 
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Since the operator iB is Hermitian, bA is real although it need not be positive. In Appendix [A| we show that for 
loa 7^ 0, the constant bA is related to the rotating- frame mode energy ea at unit amplitude by bA — £a/^>a- Finally, 
the equation of motion for the mode coefficient ca(£) is the first-order- in-time equation 



C A (t) + iLUACA(t) = T {^A 3 a cxt(*)) 
OA 



(2.37) 



Using the mode decomposition Q2.34 ), its inverse (2.35), and the equation of motion for each mode (2.37), it is 
straightforward to compute the response of the star to any externally applied acceleration a cxt (x,i). See Appendix 
[A| for justification of all the claims and formulae of this subsection. 



3. Degeneracy 

Consider next the situation where there is degeneracy, that is, where distinct modes (£ a ,uja) and (£ Bi lub) have 
identical frequencies loa = lob- We introduce the following index notation: we write the distinct eigenfrcquencies as 
to a , and we let the associated eigenvectors be £ Q k for 1 < k < n a , where n a is the degeneracy associated with u) a . 
Thus, we identify the index A with the pair of indices (a, k), and we make the identifications 



We can also rewrite sums over A as 



Consider now the matrix 



UA = LVa,k = W ) 
= to,*' 



E = EE- 

A a k=l 



(2.38) 



(2.39) 



Mab = (£ A > iB ■ Cb> + (^a + wb) (£ a > £i 



(2.40) 



Using the notation A = (a, k) and B — (b, I), we can write the matrix as Mab = M-akM- ^ n Appendix |A| we show 
that the matrix M is always block diagonal in the sense that M. a k,u = for a =/= b. In other words, for loa ^ <-^>b we 
have Mab = 0, which is a generalization of the orthogonality relation (2.32) above. 
Within a given degenerate subspace, however, we have 



(2.41) 



which is an n a x n a matrix that depen ds on the choice of basis £ a k . Since we are assuming in this section that all the 
eigenfrequencies are real, the matrix (2.41) is Hermitian and can be diagonalized by an appropri ate ch ance of basis 
£a k ~~ * ^fc so / w ithin the degenerate subspace. A choice of basis £ fe that diagonalizes the matrix ( 2.41) is the analog 
for rotating stars of orthonormal bases for non-rotating stars. For such bases, the formulae of Se c. II C 2 abo ve ar e 
valid in the degenerate as well as the non-degenerate cases, namely the mode decomposition (2.34), its inverse (2.35), 
the definition (|2.36|) of the constant bA, and the equation of motion (2.37). 



4- Reality conditions 



So far, the formalism is valid for complex Lagrangian displacements £(x, t). We now describe simplifications that 
occur in the physically relevant case of real £(x, t). 

If loa) is a solution of the quadratic eigenvalue equation ( 2. IS ), then another solution (£^, —loa) can be obtained 
by complex-conjugating the mode function and re versi ng the sign of the frequency |42[ . Under this transformation, 
the normalization constant &a flips sign, from Eq. (2.36) 0. In Appendix ^ we show that all non Jordan-chain modes 
have bA ^ 0. Hence, all the distinct non Jordan-chain modes occur in pairs (£, to) and (£*, — lo). [^] 



11 The inertial-frame frequency to + mil, where m is the azimuthal quantum numb er, als o flip s sign under this transformation. 
12 If a mode satisfies (i) £ A is purely real, and (ii) to a = 0, then from Eqs. (2.36) and (2.22) that mode also has 6a = and 
hence is a Jordan-chain mode. 
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Let us now focus attention on the set of modes with bA > 0. We introduce the following index notation: we write 
the distinct modes with bA > as (£ Q , u a ). We identify the index A with the pair of indices (a, e), where e takes on 
the values e = + and e = — , and we make the identifications £ A = £ Q e , loa = w a|e , and bA — b a , e with 



— W n 



We can also rewrite sums over A as 



5a,+ — Sc 



E = E E- 



Using these notations, the mode expansion (2.34) can be written as 



€(*) 



= E Cq .+w 



+ c Q ,_(<) 



(2.42) 
(2.43) 

(2.44) 



and its inverse (2.35) can be written as 



c a ,+ {t) = y {£ a , w„£(t) + t£(t) + *B • €(t) ) , 



and 



(2.45) 



(2.46) 



It follows t hat £ (t) will be real if and only if c Qj _(t) — c Q + (t)* for all a. In this case, writing c a = c Q + , the mode 
expansion ( 2.44 ) becomes 



kit) 



a* 



c a (ty 



t 



(2.47) 



the inverse mode expansion (2.45) becomes 

1 



c a (t) = j-{Z a , Wa€(t) + *£(*) + «B • £(t) 



(2.48) 



and the equation of motion (2.37) becomes 



c a (t) + iu a c a {t) = — (£ a oxt (i)) 



(2.49) 



Note that the number of distinct £ Q 's is the dimension of Ji. In Appendix [b] we show that the £ a 's are linearly 
independent and thus form a basis of Tt PI Hence, one comW perform a configuration-space mode expansion of the 
form 



J For modes with loa 7^ 0, we have 



b.A = oja {£ A , £ A ) H (£ a . c ■ £ a) ■ 



from Eqs. (2.18) and (2.36). Therefore positive frequency loa > would correspond to 6a > if the operator C were positive. 
However, C need not always be a positive operator. Nevertheless we suspect that for modes of nonzero frequency, positive 
frequency always corresponds to bA > 0, although we have been unable to prove this. 
14 More precisely, they would form a basis for TL if there were no Jordan-chain modes. 
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€(x,i)=5^g a (i)€ a (x), (2.50) 



instead of the phase space mode expansion _(|2.47| ). However, the coefficients q a (t) defined by Eq. (2.5C) are not related 
in any simple way to the coefficients c a (t) F~L and furthermore the equations of motion for the coefficients q a will not 
be uncoupled from one another. 

In Appendix ^ we give an alternative form of the equations of motion, valid for arbitrary complex £(x, t), which 
combines the two first order equations for c a ^ + {i) and c Q ._(t) in a single second order equation. That alternative 
form, when spec ialized to a non-rotating star and to a real basis of eigenmodes £ Q (x), coincides with the standard 



formalism ( 2.26 ) for a non-rotating star. 

Finally, in Appendix [R| we give expressions for the energy and angular momentum of perturbations in terms of the 
mode amplitudes. 



III. THE SLOW ROTATION EXPANSION 
A. Overview 



In this section we discuss the approximation method of solving for the mode functions £ A (x) and frequencies loa 
by using a perturbative expansion in the star's angular velocity SI. The main result of this section is that the leading 
order (in Q) mode functions for the inertial or hybrid modes in a slowly rotating zero-buoyancy star are orthogonal 
with respect to the inner product Q2.20D of the non-rotating star. 



There are two general classes of modes in rotating stars |14],|45j : 

• Modes for which the frequency goes to zero in the non- rotating limit. These modes can be called rotational 
since they have non-zero frequency only in the presence of rotation. For zero-buoyancy stars, these modes are 



the hybrid modes of Ref . 45 1 , also called inertial modes. For nonzero-buoyancy stars, they are purely axial. 



• Modes with a finite frequency in the non-rotating limit, such as / andp modes, and g modes in nonzero-buoyancy 
stars. 

There is a substantial literature on computing mode functions and frequencies perturbatively in powers of SI. For 
nonzero-buoyancy stars, the rotational or Rossby modes have been computed to leading order in fl by Refs. [Q, |47f| 
and p8j , and to all orders in SI by Ref. [E9| . Non-rotational modes in the nonzero-buoyancy case have been computed 



by Ref. |50|, and by Ref. |5l| beyond the leading order in SI. 

For zero-buoyancy stars, there are two classes of rotational modes. The so-called pure or classical r-modes |52j] exist 
only for I = \m\ and arc purely axial. The remaining rotational modes have both polar and axial pieces and have been 
called generalized r-modes, hybrid modes or inertial modes. They have been obtained to leading order in SI by Ref. 
[45 1 and by Ref. Q for Maclaurin spheroids, to higher order in SI by Refs. (54|] and ]5j|, to all orders in fl by Ref. 
[56 1, and in relativistic stars by Ref. |57|. The hybrid modes consist of mixtures of the r-modes of nonzero-buoyancy 
stars together with the zero- buoyancy limit of g-modes jl5) ; see also Ref. Q for hints in this direction. 

In this section we will re- derive t he eq uations obtained by Ref. [H that describe the leading order mode functions 
of the hybrid modes [Eqs (^29|)-(|33o|) below], by using a formal perturbation theory expansion of the quadratic 



eigenvalue equation ( 2.18 ). The main benefit of the analysis is that it shows that all of the rotational modes are 
orthogonal with respect to the inner product of the non-rotating star. In addition, it clarifies why, in the analysis of 
Ref. [Em, it is sufficient to consider the curl of the perturbed Euler equation and to neglect its longitudinal part. 

We also note that Ref. p(| shows that if a background star has Jordan chains of length p, and is subject to a 
perturbation of order e, then generically the leading order change in mode frequencies scales as e 1 ^ p+1 \ Now, non- 
rotating stars always have Jordan-chain modes of length 1 (see Appendix , and the leading order effects of rotation 
on the mode dynamics is linear in the star's angular velocity SI. Hence, one might suspect the existence of modes in 
rotating stars whose frequencies scale as as CI — > [jj. We show explicitly in Appendix |] that no such modes 
exist, and explain why the argument of Ref. |26|] does not apply in Appendix \ 



15 For the special case of a non-rotating star and real £(x, t), there is a simple relation between the two sets of coefficients, 
namely 2c a ,+ = q a + iq a /u a and 2c a ,- = q* a — iq^/co a . 

16 Modes whose frequencies u are proportional to at large Q have been found by Ref. |52), but for these modes cu — * luo ^ 
as f2 — » 0. 



fO 



B. Normalization convention 

The formalism of Sec. [n] is invariant under changes of normalization of the basis functions £ A ; changes of the 
normalization of the modes are compensated for by changes of the normalization of the mode coefficients CA(t). 
In general one can choose any convenient normalization convention. However, in the context of the slow rotation 
expansion, statements about the behavior of the mode functions £ A as — > are dependent on the normalization 
convention^. Throughout the rest of this paper we shall adopt the convention 

($ A ,U = MR 2 , (3.1) 

where M is the stellar mass and R the stellar radius. Then the mode coefficients are dimensionless, and order unity 
amplitudes ca ~ 1 correspond to Lagrangian displacements of order size of the star. With this convention, all of the 
mode functions £ A will have finite limits as Q —> 0. 



C. Perturbation expansion 



We denote by Ho the subspace of the Ti of Lagrangian displacements £(x) spanned by modes of the non-rotating star 
of zero frequency. In the zero-buoyancy case, the space Tto consists of those perturbations £(x) for which 5p = — V-(p£) 
vanishes, and in the nonzero-buoyancy case it consists of purely axial vectors (under reasonable assumptions). These 
characterizations of the space of zero- frequency modes have been proved for fully relativistic stars by Ref . |57j| . In 
Appendix ||| we give Newtonian versions of the proofs of Ref. 1 57 , and also derive some other properties of the space 
Ho- 

The operators B and C can be expanded as power series in the angular velocity Q as [£] 

B = f>B (1) (3.2) 

and 

c = c<°> + n 2 c^ + 0(n A ). (3.3) 



Also the inner product ( 2.2C ) depends on f2 through the background density p and can be expanded as 

<^O = (£,O + tt 2 (£,O 2 + o(ft 4 ). ( 3 - 4 ) 

Since the operators C and iB are Hermitian with respect to ( , ), it follows that and iB^ (although not C' 2 ') 
are Hermitian with respect to ( , ) . 

The operator C' ' governs the modes of the non-rotating star. Let u>^\ for A = 1, 2, 3 . . ., be a complete set of 
solutions of the eigenvalue equation for the spherical star: 

c^-L = «TL- (3-5) 



As in Sec. [I C 4 above we can take A — (a, e) with e = + or e = — , and with £ a + = £ a _ = £ a , and = — = 

LJa^ > 0. Since C^ ^ is a Hermitian operator, the basis {^ Q } is a complete, orthonormal basis of TL. The basis £ a 
can be chosen arbitrarily within each degenerate subspace of the operator C'- -'. For example, the basis can be chosen 
arbitrarily within the subspace H.Q- Below we will introduce a different, preferred basis of eigenvectors of C^ -* 
that are the — > limit of modes of the rotating star. 

Consider now trying to solve for the one parameter family or families of modes £ a (fi),a; a (Q) of the rotating star 

for which u> a (Q,) — * uja^ as O — ► PL As explained above, there are two types of such modes, rotational modes for 



17 For example, if one chooses to normalize an r mode in such a way that an amplitude ca of order unity corresponds to a 
mode energy of order the stellar binding energy, then £ A oc as Q — > 0. 

18 More precisely, the expansion parameter is the dimensionless quantity f2^/ _R 3 /(GA/), where M and R are the stellar mass 
and radius. 

19 It suffices to consider the case e = + and to seek a one-parameter family of modes whose limiting frequency is since 
a on e para meter family whose limiting frequency is can then be obtained using the transformation (£, oj) — > (4* 
Sec. 



1IC4 above. 
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which LOa^ — 0, and non-rotational modes for which lu^ ^ 0. For both types of modes, we make the ansatz that the 
frequency and mode function can be expanded as 

cj a (Q) = 4 0) + r*4 1} + ft 2 + o{n 3 ), 

and 



(3.6) 



(3.7) 



If we substitute the expansion s (j3li| )-( |3.7D and also the expansions ( |3.2| )~( |3~3| ) for the operators B and C into the 
quadratic eigenvalue equation (2.18) we get a series of equations. First, at order O(f2 ) we get 



. ro) 2 + c (o) 



f(0) 



0, 



which is the usual eigenvalue equation for spherical stars. At order 0(f2) we get 



_ W (0)2 + C (0) 



Finally, at order 0(f2 2 ), we get the equation 



(0) 



_ 6L ,(0)2 +C (0) 



(0) 



2 W ! 1 >+ i B( 1 ' 



2 W W+ i B( 1 



-2w(%«-a;« a -ia;WB« + C« 



• d 0) = o. 
d 0) - o. 



(3.8) 



(3.9) 



(3.10) 



We next express these perturbation equations in the basis of modes of the non-rotating star. We expand the mode 
functions an d Ci 2 ' m terms of the eigenbasis {$, a } of the non-rotating star as 



;(0) _ 



E C af3 £e 

P 



(3.H) 



(3.12) 



and 



e(2) = 

Sft 



Vc (2) 

/3 



(3.13) 



"a/3 = *W because of the possibility of degeneracies. Just as in degenerate 
perturbation theory in quantum mechanics, the modes of the rotating star define a preferred basis £^ of each 
degenerate subspace. Since we do not know this basis before we solve for the modes of the rotating star, we must 
start our computation with some arbitrary choice of basis £ a in each degenerate subspace. 

Next, we substitute the expansion (pl| ) of ^ into the zeroth-order eigenvalue equation ( |3.8| ). This gives, using 
Eq. (ph, 



E^ (0) 



a/3 



4 0)2 -4 0)2 



€/3 = o. 



(3.14) 



It follows that cl°2 = when oj a °^ ■ In other words, the matrix cfl — is block diagonal, with each block 



corresponding to a degenerate subspace. Hence the zeroth order mode function can be written as, from Eq. (3.11) 



t(0) = 



E 



JO) i 
c a/3 <5/3' 



(3.15) 



and thus lies within the degenerate subspace corresponding to the eigenvalue tOa \ In particular, for the modes with 
(iif = 0, the zeroth order mode function £^ lies inside the space Ho- 

The analysis now divides into two cases for the two different types of mode discussed above. 
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D. Modes with finite frequency in the non-rotating limit 



If we substitute the expansions (3.11) and (3.12) for £^ and into the first-order eigenvalue equation (|3~9|), and 
then take the zeroth order inner product ( , ) with we get 



(0)n2 



i 1 ) , O u (0) (1) (0) 



,(0) 



V r(°) - 



(3.16) 



where 



If we specialize Eq. ( |3.16| ) to values of for which wl ^ = w£° , and divide across by we get 

2 W (1 »f (0) + iVB (1) c (0) - 



(3.17) 



(3.18) 



which is a standard eigenvalue equation. Thus, the zeroth order modes are given by diagonalizing the operator 



within the degenerate subspace of modes (3 with oA "* = u>a' ■ The eigenvalues of give the first order changes uj, 
to the eigenfrequency, just as in degenerate perturbation theory in quantum mechanics. 



,(0) 



»(*) 



01 
(1) 



Suppose now we switch to the basis that diagonalizes B^l within each degenerate subspace, so that 
Then, applying Eq. (3.16) for values of {3 with uj^ ^ oj^ gives for the first order change in the mode functions 



JO 



lap- 



■ (0) 



B 



(i) 



,(0)2 (0)2 Pa' 

U I3 



(3.19) 



Using the expansions (3.7), ( |3.11 ), and (3.12) together with Eq. (3. IE) shows that the inner product of two mode 



functions with LUa^ ^ UJ^ 1 is 



■in 



(0) 



B 



(i) 



(0) a/3 




o(n 2 ). 



(3.20) 



This explicitly demonstrates that modes in a rotating star are not orthogonal in general, since we can find pairs of 
modes a, /? of the non-rotating star for which the matrix element B^ is non vanishing. 



E. Modes with vanishing frequency in the non-rotating limit 

Things work somewhat differently for the second class of modes. If we substitute uj^ = into Eq. (3.9) we get 

C (0) . = 0. (3.21) 

Thus, the first order equation ( |3.9| ) does not determine \ except to dictate that £^ lie in the space Ho of zero- 
frequency modes of the non-rotating star. In addition, the ambiguity in the choice of basis of £^ ' of TCq is not resolved 
at this linear order. 

(0) 



Hence we must use the second order equation (3.10). For u>a = 0, this equation reduces to 



C (0) -d 2) + -ai^-^BW+cP) 



(o) _ 



0. 



(3.22) 



If we substitute the expansions ( 3.11 ) and ( 3.1 3| ) into Eq. ( |3.22 ) and take the zeroth order inner product ( , ) with 
C/3. w e get 
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1 



(2) 
01 



c ( °) = 



where 



(7 (2) = (i ■ £ 



(3.23) 



(3.24) 



We now specialize Eq. (3.23) to values of /? for which a;^ "* = 0, i.e., project both sides of the equation into Ha 
This gives 



4 1)a c«3+E[-^ 1)fl S 



) , r (2) 



.(0) _ 

~ory 



(3.25) 



which is a quadratic eigenvalue equation that determines both Wa and c^i. Equation (3.25) can be rewritten as 



Po 



Po • = 0, 



(3.26) 



where Po is the orthogonal projection operator (with respect to ( , ) ) that projects into H.q (see appendix]^). 

It is at this point in the analysis that the difference between stars with zero and non-zero buoyancy enters. We 
treat the two cases separately. 



1. Zero-buoyancy stars 
When A = 0, the term P • ■ P • ^£ 0) in Eq. fl3.26| ) vanishes. To see this, let 

K =c.e=c.p .e 



(3.27) 



It follows from the formula ( |2.13 ) for the operator C that k(x) must be a pure gradient, so that Vxk = 0. However, 
Eq. (E6) of Appendix ^ then implies that Po ■ ac = 0, and expanding this last equation to second order in fl gives the 
desired result. 

Equation ( ^26| ) now simplifies to 



p -;bW-p 



;(0) 



(i) t(o) 



(3.28) 



which is a standard eigenvalue equation in Tio- Thus, the zeroth order mode functions £^ ' f° r the hybrid modes 
are simply the eigenvectors in Hq of the Hermitian operator iPo ■ B^ • Po- It follows that these zeroth order mode 
functions are orthogonal to each other with respect the inner product ( , ) of the non-rotating star. They are 

also orthogonal to the zeroth order mode functions of the modes with ^ 0, since those mode functions £^ lie 
completely inside degenerate subspaces that are orthogonal to TCq, by Eq. ( |3.15 ). 

We now rewrite the equations defining the hybrid modes £^(x) in a more accessible notation. First, from Eq. 
( |3 .15 ) we know that lies in Hq, which from Appendix |] is equivalent to 



(Oh 



0. 



Second, Eq. ( 3.28 ) can be rewritten using the definition ( |2.8|) of B and Eq. (E6) as 



V x 



(0) 



(3.29) 



(3.30) 



Equations ( |3.29 ) and ( 3.30| ) are the equations that are solved to obtain the hybrid modes in Ref. ]45j] . 



The component of Eq. (3.23) orthogonal to TCo is not needed to determine the leading order quantities ' and wl, 



It 



determines the higher order correction to the mode function 
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Finally we note that going from Eq. (3.26) to Eq. ( 3.28Q entailed dividing by u>a ■ There is another solution to Eq 
( ^26|) with 



,(i) - 



0. Examining the perturbation expansion to higher order for this solution shows that u, 



(2) _ 



also. We suspect that this solution is a zero-frequency solution to all orders in Q. These modes are not Jordan-chain 
modes, but they are pure gauge in the sense of Ref. [|IJ, at least to linear order in the mode amplitude. We note 
that it is impractical to include these modes (and also the similar modes of Appendix in nonlinear mode evolution 
calculations since the perturbation expansion in mode amplitude for these modes typically breaks down on a timescale 

~i/a 



2. Stars with buoyancy 



For A / (stars for which there is a nonzero buoyancy force), the space Ho is the space of purely axial vectors 
(see Appendix ^). In this case the quadratic eigenvalue equation ( 3.26 ) for the mode functions can be written as 



Po 



,(i)2e(o) 



(0) 



p(2) AO) 



0, 



(3.31) 



(21 (2) A 

where Po • v is now the axial part of v, and C„ and Cj, are the second order in ft pieces of the operators ( 2.13| ) 
and (2.15). The term • £^ drops out since it is a gradient and thus has no axial part, and the term c[ 2 ' • £^ 
drops out since it depends only on V • £^ ' and • e r which vanishes as |i 0) is purely axial. Using Eq. flEll| ), the 
equation now reduces to 



r ■ V x 



(0) 



(3.32) 



In Ref. E^j it is shown that Eq. (3.32) determines the frequency uj^ but not the eigenfunction £ffl , and that one 
must go to higher order in the Vl expansion in order to determine the eigenfunction; see also Refs. J46| and [Q. 



IV. SECOND ORDER LAGRANGIAN PERTURBATION THEORY IN ROTATING STARS 



In this section we develop a second order Lagrangian perturbation theory for rotating stars, extending previous 
analyses of non-rotating stars (see, e.g., Refs. | p9|j60| , p8[ ) . The m ain n ew result is that the expression for the three- 
mode coupling coefficient in terms of the mode functions [Eq. ( 4.20 ) below] in the same for rotating stars as for 
non-rotating stars. Although we restrict attention in this paper to second order perturbation theory, we lay the 
foundations in Appendix |j| for computing mode coupling coefficients at third and higher orders. 



A. Second order equation of motion 



To describe the fluid perturbation we use the rotating-frame Lagrangian displacement £(x, t). The definition of this 
quantity is the same in a nonlinear context as in linear theory: a fluid clement at rotating-frame location x in the 
unperturbed star is moved to rotating-frame location x + £(x, t) by the perturbation. We assume that the pressure 
p in the background star is a function p(p, ft) of the the background mass density p and of a vector p of other fluid 
variables such as entropy or composition, which can vary with position. We also assume that the fluid variables p of 
each fluid element are preserved by the perturbation. 

Under these assumptions, the equation of motion for £ is derived in Appendix | to second order in £ via direct 
perturbation of the Euler equation, and in Appendix [J] to any order in £ from a variational principle. The result has 
the form 

£ B i C i a[£], (4.1) 

where the left hand side contains all the linear terms [cf. Eq. ( p. 7] ) above] and the right hand side a[£] is a nonlinear 
acceleration. That acceleration can be written as a = ap + &g, where &p and are the contributions from pressure 
gradients and gravity; the corresponding force densities are f = pa, fp = pap, and {q = P&g- 
We expand the acceleration as 

= a (2) K, + 0{e) = 4 } K, £] + 4 ] K. £] + 0(H 3 ), (4.2) 
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where a^ 2 ) [£,£'], ap 2 -* [£,£'], and [£,£'] are symmetric bilinear functions of their arguments. In order to give 



the explicit results for the second order pieces a P [£,£] and a£j [£,£] of the pressure gradient and gravitational 
accelerations, we need to introduce some notations. The first and second order Euleri an p erturbations 5^ <fi and 5^(j) 
to the Newtonian potential are defined as functionals of £ by the equations [see Eq. (130) below] 



and 



J_V 2 <5( 2 )0=iv i V j ( P ee i )- 



We define the tensors 1 ■ , S l ■ and x 1 by 



(-)' 



fe4 ;/■ 



(4.3) 



(4.4) 



(4.5) 
(4.6) 



and 



(4.7) 



Here semicolons denote covariant derivatives, for example ■ = Vj£\ Fin ally, we denote by Li the generalized 
adiabatic index governing the perturbations, defined b y Eq . ([16) or Eq. ( [J24| ). 
The second order gravitational acceleration is [Eqs. ( |l37| ) and ( [J32| ) below] 



(4.8) 



The second order pressure term is [Eqs. (137), (J23) and (JJ28)) below] 



P (ri-i)e%+pH J 



(4.9) 



where 



-p6 



(ri 



D 2 



9T! 



91np 



:P(ri 



IB, 



(4.10) 



and O and S are the traces of J , and S J ^ . More general expressions that are valid to all orders in £ are given in 
Eqs. (|J2|) below. 

A key feature of the equation of motion fl4.l| ) is that there is no dependence on the angular velocity fl of the rotating 
frame in the nonlinear terms on the right hand side. The dependence on fi occurs only in the linear terms on the 
left hand side. This can be understood most easily from the variational principle that underlies the hydrodynamical 
equations. As explained in Appendix [j], the only term in the action that is not explicitly invariant under a change of 
rotational frame is the kinetic energy term, and that kinetic energy term is quadratic in the Lagrangian displacement 
£ and its time derivative £. 

Finally, we note that in Lagrangian perturbation theory only pressure and gravity can contribute to the nonlinear 
force. In spite of the fact that the nonlinear advection term v • Vv appears in the Euler equation, no nonlinear 
advection terms are present in the Lagrangian framework, since the kinetic energy is explicitly quadratic in 



B. Coupled equations of motion for the mode coefficients 



We now use the equation of motion given by Eqs. (4.1), (p~^), (4.8) and ( p~9| ) to derive the coupled, nonlinear 
ordinary differential equations satisfied by the mode expansion coefficients ca$)- The basic equation (2.7) of Sec. || 
can be applied with the externally applied accelerat ion a ext (x, t) replaced by the nonlinear acceleration a[£]. Therefore 
we can directly carry over the equation of motion (2.37), which gives 



c A (t) + iiu A c A (t) = (£ A , a( 2 > [£, £] ) + O(e)- 



(4.11) 
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Next, we take the complex conjugate of the phase space mode expansion ( [2.34 ) and use the fact that £(x, t) is real to 
obtain 



£(x,*) = £c4i)*uxr 



(4.12) 



Substituting this expansion into the right hand side of Eq. ( 4.11 ) gives 

t A {t) + iLu A c A {t) = j-^2 k* abc c B (t)*c c (ty + 0(c 3 ). 



B,C 



where the three-mode coupling coefficient is 



^ABC = {C A ,^bAc. 



which is completely symmetric in the indices A, B and C . For oja ^ 0, this equation can also be written as 

c A (t) + iu A c A {t) = iu A V c B {t)*c c (ty + <3(c 3 ), 



B,C 



£A 



(4.13) 



(4.14) 



(4.15) 



where we have used Eq. (K22). This form of the equation has a simple physical significance, since kabc/^a is the 
ratio of the nonlinear intera ction energy at unit amplitude to the energy of the mode at unit amplitude. 

The equation of motion (4.13) is valid for arbitrary complex £(x, t). However, for the physically relevant case of 
real £(x, t), the coefficients CA(t) will not all be independent, as explained in Sec. [IC4 above. The modes will occur 
in complex conjugate pairs (£,w) and (£*,— lo), and it is convenient to use the set of modes £ Q defined in Sec. II C 4 
consis ting of one mode from each complex-conjugate pair. For this set of modes, the mode expansion is given by Eq. 



€(t) = £M*)€« + q«(*)*£] 



(4.16) 



and the inverse mo de ex pansion by Eq. (2.48). For this e xpansion, the variables c a (t) are al l ind ependent. If we now 
use the expansion ( 4.16] ) in the equation of motion (2.49) we find as a replacement for Eq. (4.13) the equation 



L a: 1 a 



I 

b a 



where 



/3,7 



(4.17) 



(4.18) 



which is symmetric in a, (3, and 7. In Eq. (4.17) we have used a notational convention for the coupling coefficients k 
where a bar over an index means that the corresponding mode function is to be complex conjugated in the expression 
for the coupling coefficient. For example, we have 



(4.19) 



Equation (4.17) is the final equation of motion. To model the saturation of r-modes, it will have to be supplemented 
by viscous damping and gravitational radiation reaction terms. We now turn to evaluating the coupling coefficients 
kabc to insert into this equation, with the understanding that each capital Roman index can b e either a unbarred 
or a barred Greek lower case index, for example A — a or A — a. In the notation of Sec. II C 4, (a, +) corresponds 
to a and (a, — ) corresponds to a. 



C. Explicit expression for the three-mode coupling coefficient 

We can obtain an explicit expression f or t he cou pling coef ficient kabc i n terms of the mode functions (x) , £ B (x) 
and ^(x) by substituting the formulae (4.2), (4.8) and (4.9) for the nonlinear acceleration into the definition (4.14) of 
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kabc, using the definition (2.20) of the inner product, and performing several integrations by parts. The computation 



can be simplified by (i) computing the functional (£* , a*- 2 ' of (ii) equating this functional to £) where 

,£ ) i s a symmetric trilinear functional of its arguments, in order to determine the functional £ , £ ), and 
(iii) evaluating k(£ai£b>€c) to obtain the coupling coefficient kabc- For the gravitational terms, the dom ain of 
the spatial integral can be taken to be all of space, if one includes the ^-function contributions in Eq. ( |I30|) . The 
gravitational terms can then be integrated by parts, and the boundary terms at infinity which are generated vanish. 
The domain of integration for the pressure terms can be taken to be the interior of the star; the boundary terms 
generated by the integration by parts vanish since p = on the stellar surface. The final result is 



kabc = \ I 



(Ti - 1) 



(T a - 1) (E^s V • i c + E BC V • £ A + ScmV • £ B ) 
dT 1 



d\np 



V • £ A V ■ £ B V • £ c + XABC + XACB 



~ 2 J d3 *P l&tB^faiJ +Ck^ (1) ^A ; « +&& {1) <t>B;i3 



(4.20) 



The notations here are as follows. The quantity 8^4>a is defined by Eq. (4.2) with £ on the right hand side replaced 
by $,a- We define the functions 



and 



The tensors 1 j , 5 l j and \ l j which we defined previously in Eq. 
at £ = £ = Finally we define the traces 



and 



~AB=51 Z l M A ,U 



XABC = &iX l j[^A^B,ic\- 



(4.21) 
(4.22) 

(4.23) 

are obtained from these expressions evaluated 

(4.24) 

(4.25) 



The domain of integration for the gravitational term in Eq. (4.2C) is all of space. The domain of integration for the 
pressure terms can be taken to be either the interior of the star, or all of space. 

Note that the expression ( 4.20 ) for the three mode coupling coefficient is valid for arbitrary angular velocities VI 
of the background star. The expression has n o explicit dependence on Q, however, and therefore is also valid for 
non-rotating stars. In fact the expression ( 4.20 ) agrees with previous results |^,^,^8| for non-rotating stars, except 
for an overall factor of —3 due to differing conventions for the form of the equation of motion. Some of the previous 
analyses omitted the pressure force term containing the factor dTi/dliip, and some used the Cowling approximation 
and thus omitted the gravitational terms^J. 

For the special case of constant Ti, it is possible to obtain an alternative form of the pressure terms in the coupling 
coefficient by (i) taking the domain of integration in both the gravitational and pressure terms to be all of space, and 
(ii) performing further integrations by parts pl| . The result can be expressed as 



K ABC = k(£a>£bi£c) 



(4.26) 



where 



21 Ref. [[jo| makes a slight error in calculating the gravitational terms. However, since they use the Cowling approximation in 
which SS — 0, this error does not affect their results. 
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«(€,€,€) = -§ / d 3 x 



eez k P;iik + 3(v • oeo^ - > . = i , + imv • o a 



+ 3f + e(v ■ seisph + 3et 3 pm^ + ree fe / 



(4.27) 



and (5p = —T\p (V ■£)—£■ Vp is the Eulerian pressure perturbation. The domain of integration here is all of space, 
and therefore one should include the 8- functions at the stellar surface in the factors p-jjk, P-ij and possibly <fi-,ijk- 
However, one can show that all of the contributions from the S functions cancel in general, and therefore one can take 
the domain of integration in the expression (4.27) to be the interior of the star. 



D. Application to Rossby Modes 



Although the expression for the three-mode coupling coefficient is the same for rotating stars as for non-rotating 
stars, it is much more difficult to evaluate the coupling coefficient for rotating stars, for two reasons. First, including 
centrifugal flattening of the background pressure and density profiles makes the integrals much more difficult to 
perform. Second, the mode functions £ A (x) for rotating stars typically consist of a sum of terms with different values 
of the quantum number I, which complicates the calculation. 

We shall avoid the first of these difficulties by using the slow rotation expansion. First, we will compute the 
mode functions to leading order in Q, as explained in Sec. Ill above. Second, in evaluating the coupling coefficient 



expression (4.20), we will use the background density, pressure and potential profiles of a spherical star, i.e., we will 
use the — ^ limits of the variables of the background star. This procedure will give the zeroth order piece of the 
coupling coefficient kabc in an expansion in powers of f2. In cases where this zeroth order piece vanishes, the coupling 
coefficient might be nonzero at higher orders in il. 

Before discussing in more detail how to compute the coefficients, we describe some useful selection rules. 



V. SELECTION RULES FOR THE COUPLING COEFFICIENTS 



A. Terminology 



In this section we derive a number of selection rules that apply to the three- mode coupling coefficient. We start by 
reviewing some terminology to describe modes. Any vector field £(x) can be expanded in vector spherical harmonics 

as 



! in 



Ai m (r)Y lm (6, tp)et + Bi m (r)VY im (9 , (p) + C hn (r)rxWY lm (9, tp) 



(5.1) 



where (r, 9, p) are spherical polar coordinates and and is the associated orthonormal basis. Modes for which 

Aim = Bi m — for all I, m are called axial modes. We define h to be the parity map [x, y, z) — > (—x, —y, —z) and / 
be the map (x, y, z) — * {x,y,~z) or (r,9,ip) — > (r, ir ~ 9,tp) which we call z-parity. Rotating stars, both Newtonian 
and relativistic, are invariant under these maps, and therefore we can always choose to use a basis of modes for which 
all the modes have definite parity and z-parity transformation properties; see Refs. p5| , |62| . Specifically, we define the 
pullback /»£ of £ under / by |63| 



/.€ = /* [T (r, 0, p)e f + e(r, 9, <p)e § + ^(r, 9, V )e^ 

= C (r, tt-9, <p)et - 7T-6, p)e § + &(r, tt - 0, <p)e . 
Then we can always find a basis for which all modes are either z-parity even, satisfying 

M = €, 

or z-parity odd, satisfying 

/*£ = -e 



(5.2) 
(5.3) 
(5.4) 
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Note that f*Yt m = (— 1) m Yi m , A similar statement is true for the regular parity map h; we can take all modes to 
be either parity even or parity odd. Finally we can take all modes to have definite values of m, so that if z is the 
mapping (r, 9, <p) — > (r, 0,<p + Aip), then 

z»£ = cxp[imAip}$,. (5-5) 
Recall also that all mode s can be classified as either rotational (u> — > as — > 0) or regular (uj — > finite as Q — > 0), 



as discussed in Sec. [II A above 



B. Classes of modes 



We next review the various different classes of modes to which we shall apply the selection rules. In zero-buoyancy 
stars, we distinguish three different classes of modes: 

• The rotational modes which are axial to zeroth order in Q. These modes are characterized by single values of 
/ and of m. If one is working with all the modes then one can have m = I or m = —I. If one is working 
only with the set of modes £ Q with positive rotating-frame frequency u) (the remaining modes being complex 
conjugates of these), then only m = — I is allowed. These modes are all z-parity odd. We shall call these modes 
the pure r-modes. 

• The rotational modes whose O — ^ limits are not axial. These have fixed m but do not have a fixed value of I. 
To zeroth order in f2 they can be expanded as |45| 

00 

CM = E 

3=0 



W j+m+1 (r) 



Y™ m+1 (r)e f + V j+m+1 (r)VY% m+1 



j-\-m 



(5.6) 



in the z-parity odd case, and as 



~W j+m (r) 



Y™ m e f + V 1+m (r)VY, 



U J+m +i(r) 



LYl 



(5.7) 



in the z— parity even case, where L = — i(r x V) and j runs over 0,2,4, . . .. Following Ref. |57| we shall call 
these the rotational hybrid modes, although they have also been called inertial modes or generalized r-modes. 

• The regular modes, which can be z-parity even or z-parity odd, and consist of / and p modes. 
Similarly in stars with buoyancy we shall distinguish the following classes of modes: 

• The rotational modes. All of these modes are axial to zeroth order in f2 (Appendix They are also character- 
ized by single values of I and m, but any values of I, m with I > \m\ are allowed. They have z-parity (— iy+ m + 1 
and parity (— 1)' +1 . We shall call these the r-modes. [Unlike the pure r-modes above, their radial eigenfunctions 
are not polynomial in r.] 

• The regular modes, which may be z-parity even or z-parity odd, and consist of /, p and g modes. 



C. Selection rules 



We can write the three mode coupling coefficient as a symmetric function of the three mode functions £^4, £ c : 

habc = k(£ a , £ b ,€c)- ( 5 - 8 ) 

Now apply the pullback /» to this equation. Since the left hand side is a pure number, it is invariant. On the right 
hand side, we can pull the /* operator inside the function k since pullback commutes with geometrical operations like 
taking covariant derivatives. Thus 

kabc = «(/*£a> /*£bi /*£c)- ( 5 - 9 ) 
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If we denote the z-parity of mode A by ea = ±1, so that J*£a = € a£a> then it follows from Eqs. (5.8) and ( |5.9| ) that 



[I - e A eBec] kabc = 0. (5.10) 

Therefore we get the selection rule 

z — parity selection rule : odd number of z-parity odd modes =>■ kabc — 0. (5-11) 

For example, if all three modes are z-parity odd, then the coupling coefficient vanishes. An identical argument gives 
the selection rule for regular parity 

parity selection rule : odd number of parity odd modes => kabc — 0. (5-12) 

Fina lly, applying Eq. (|5.10|) with / replaced by the mapping <p — > tp + Aip and with 6a replaced by exp[im,tA</3] [cf. 
Eq. ( (5.5| ) above] yields 

m selection rule : nriA + tub + mc ^ =>■ k^bc = 0. (5.13) 

These selection rules are valid in rapidly rotating stars as well as in the slow rotation limit, and also retain their 
validity for relativistic star s. The ar gum ents can also be generalized to four-mode coupling coefficients kabcd- 

The selection rules ( pUl] ), ( fjf ), (£t|) are based on symmetries of the background star and thus apply to all orders 
in d. Therefore, when they are applicable, they restrict not only the piece of kabc which is zeroth order in d (which 
is all we compute in this paper), but also all correction terms that are higher order in d. By contrast, in Appendix 
p| we derive two selection rules which are valid to zeroth order in d only. First, we show the coupling coefficient for 
the coupling of three axial modes is 0(d). Thus we have 

Axial selection rule 1 : $, B , £ c = axial + 0(d) kabc — 0(d). (5-14) 

We also prove in Appendix ^ the following second selection rule. If (i) two of the modes are axial, (ii) the third has 
vanishing Eulerian density perturbation, and (iii) the background star and perturbations obey the same one-parameter 
equation of state (i.e., the zero-buoyancy or barotropic case), then kabc — 0(0). We can write this selection rule as 

Axial selection rule 2 : bp a = 0(d), £b>£c = axial + 0(d), zero buoyancy =4> kabc — 0(d). (5.15) 



Note that axial modes automatically have vanishing density perturbation to zeroth order in d, from Eq. (2__3) and 



using V - £ = C = and the fact that the background density is a function of r only to leading order in d. The 



selection rules ( 5.14 ) and ( 5.15 ) might fail for relativistic stars. 

In incompressible stars, we prove in Appendix |G] that the coupling between any three rotational modes vanishes to 
zeroth order in 0. A precise statement of this selection rule is 

Incompressible selection rule : SpA, Spa, 5ps, $Pb, Spc, Spc — 0(1/Ti) + 0(d) kabc = 0(1/Ti) + 0(d), 

(5.16) 

where Ti is the adiabatic index of the perturbations. Note that the regular modes (/-modes) of an incompressible star 
have 5p — > const in the limit Ti — > oo, and do not satisfy the hypothesis of this selection rule, whereas the rotational 
modes have exactly vanishing Sp and dp (to zeroth order in d) and so do satisfy the hypothesis. 

The d dependence of the coupling coefficients kabc arises in two ways: (i) through the dependence of the expression 



(4.20) on the background star, for which the leading order corrections arise at 0(d 2 ), and (ii) through the dependence 
of the mode functions themselves on d, for which the leading order corrections are 0(d 2 ) for pure r-modes |54| , [5"5"f , 
but can be 0(d) for regular modes. Therefore, when the zeroth order coupling coefficient vanishes, the leading order 
corrections may be 0(d) or 0(d 2 ), depending on the modes involved. 



D. Applications of selection rules 



We now discuss some applications of these rules. Consider first the case of zero-buoyancy stars. The most interesting 
coupling coefficients are those in which an unstable pure r-mode appears twice, since if the amplitudes of all the other 
modes are zero initially, it is via these coupling coefficients that energy can leak out of the pure r mode. In particular 
we are interested in the most unstable mode, the I = m = 2 mod e. If £ A = £ B is an unstable pure r-mode, then 
the third mode £ c must be a z-parity even mode from rule (5.11). In particular, the third mode cannot be a pure 
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r-mode as pure r modes have odd z parity. Coupling coefficients between any three pure r-modes vanish. Nonzero 
coefficients can be obtained by taking the third mode to be a z-parity even regular m ode. However, if the third mod e 



is taken to be a ^-parity even rotational hybrid, then the second axial selection rule (5.15) applies due to Eq. (3.29), 
and implies that the coupling coefficient vanishes to zeroth order in £1. It may be nonvanishing to higher order in f2 
p4| . If one considers one pure r-mode and two rotational hybrid modes, the corresponding coupling vanishes in an 
incompressible star by the incompressible selection rule, but may be nonvanishing at 0(1/Ti). 

In Table Q we summarize which classes of modes can and cannot couple to zeroth order in Q in in zero-buoyancy 
stars. Stars indicate non-zero coupling coefficients. Since kabc is symmetric, Table [j] lists all possible types of coupling 
coefficient except one, the coupling coefficient between one pure r-mode, one regular mode and one hybrid rotational 
mode; this coupling coefficient is non-zero. The table shows the kinds of coupling coefficients that can be non-zero 
when z— parity restrictions are met. For example, the table indicates that three hybrid rotational modes can couple 
together, but this can only happen when the coupling involves an even number of z— parity odd modes. 

For nonzero-buoyancy stars, consider the coupling coefficient between three r-modes. If the three modes have 



quantum numbers Ia, rriA, Ib, ms and lc, m c, then since each mode has z-parity (— l)' +m+ , Eq. (5.1C) gives 



[l - (_l)H»>A+m B +m +U+! B +l c ] Kabc = q ( 5 17 ) 



Invoking the m selection rule (5.13) implies that Ia + Ib + lc must be odd for kabc to be non-zero. If all three modes 
are of the \m\ — I variety, this contradicts tua + tob + mc = 0, so the coupling of three I = \m\ modes is zero. If some 
of the three modes do not satisfy I = \m\, then the ^-parity and parity selection rules do not apply. However, the 



axial selection rule (5.14) does apply and therefore the coupling coefficient vanishes to zeroth order in f2. To higher 
order in Q however, these coupling coefficients with I a + Ib + lc °dd and ttia + rriB + mc = need not vanish. In 
fact, Morsink [Q has calculated coupling coefficients of this type to 0(Q 2 ) and has found them to be nonzero. 

VI. COMPUTATIONAL METHOD FOR COUPLING COEFFICIENTS 

In this section we describe an efficient computational method for computing coupling coefficients in rotating stars. 
Specific coupling coefficients relevant to r-mode saturation will be computed numerically in a subsequent paper. 

A. Overview 



Direct calculation of the expression ( 4.2C ) for the coupling coefficients is algebraically intensive even in the case 
where all three modes have no sums over I. If one proceeds in spherical coordinates to directly calculate all the 
covariant derivatives needed one encounters many terms that scale as 1/ sin™ 8 with n > 2. These terms arc divergent 
as 9 — ► or tt. Since the full coupling coefficient is obviously finite at the poles (9 = 0,7r), these terms must cancel 
with other terms of the same order in 1/ sin#. However, the amount of algebra necessary to cancel these terms by 
hand is unmanageable. Calculation of the coupling coefficients in Cartesian coordinates is also daunting. Even though 
all covariant derivatives are highly simplified in Cartesian coordinates, the total number of terms is greater than in 
spherical coordinates and one must still integrate over the sphere in the end, which is awkward. Cartesian coordinates 
also obscure selection rules that might lead to insight into which modes couple and which don't. Methods based on 



integration by parts simplify the algebra for non-rotating stars |28 6CJ, but, for rotating stars, modes acquire toroidal 
pieces that severely complicate the angular integrations even after integration by parts. 

To circumvent these difficulties, we have found that using the spin-weighted spherical harmonics of Ref. p^j] as 
angular basis functions is ideal. The use of these basis functions allows one to immediately reduce the angular 
integrals to Wigner 3 — j symbols. It also reduces dra matica lly the number of terms that must be calculated. For 



example, for modes that have only one value of I each, ( 4.20 ) requires explicit calculation of ~ 100 terms, including 



many that require computing angular derivatives of the mode functions. When modes involve sums over n values of I, 
the number of terms to calculate scales as n 3 . Once this is all done, one must then integrate the resulting complicated 
expression over angles. However, in the spin-weighted formalism, the calculation involves just 6n 3 terms, none of 
which requires computing angular derivatives. Moreover, those 6n 3 terms are completely spherically symmetric, all 
angular dependence being already integrated out in terms of Wigner 3 — j symbols. The coupling coefficients are then 
built from these spherically symmetric terms, their corresponding 3 — j symbols, and their complex conjugates and 
symmetrizations. The method is easy to automate, including the various sums over the different Vs of the modes, 
leaving only 6 terms to actually code up. 



This section is organized as follows. In Sec. VI B we define the spin-weighted basis and express the covariant 



derivatives needed for the coupling coefficients in terms of components on that basis. Some details are relegated to 
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Appendix ||. Section VIC derives the explicit expression (6.24) for the coupling coefficient in terms of spherically 
symmetric functions and Wigner 3 — j symbols. 

The method described here is specialized to slowly rotating stars and is valid only to zeroth order in f2, since we 
take the background density and pressure to be functions of r only. However, it is straightforward to extend the 
method to higher orders in Q. For example, in calculating coupling coefficients to second order in Q, contributions 
of order (9(S1 2 ) will arise both from corrections to the background stellar model (centrifugal flattening), and from 
corrections to the modes. Including the additional angular dependence of the background density and pressure would 
require computing angular integrals of products of four spin-weighted spherical harmonics instead of three. 



B. An appropriate choice of basis 



Consider a given mode function £(x). Since any mode has a fixed value of to, we can assume an expansion of the 
mode function on a basis of vector spherical harmonics of the form 



l — \m\ 



Wi m (r) Uim(r) TV 

Y lrn e f + Vi m (r) vYim LY h 

r r 



(6.1) 



where L = -i(r x V). We define a set of spin-weighted basis vectors consisting of two complex null 3- vectors, m and 
m given by 



and 



m = 71 ^ + 



and a vector orthogonal to those two, 

1 = e,.. 

The components of the metric on this non-coordinate basis are 

9 U = 911 = 1, 

y tJmm J-j 

with all other components zero. The components of the mode function £ on this basis are given by 

6 = #o = I • € 

£m = i?i = m • £ 
Cm =•#_!= m • £, 



(6.2) 



(6.3) 



(6.4) 



(6.5) 
(6.6) 



(6.7) 



where the subscript 0, 1 or —1 on the scalar function z9(r, 0, tp) represents the spin weight^ of the component. Note 
that 



6 = = 

6n = C = 01 

and that the mode function can be reconstructed as 

£ = C'l + Cm + Cm. 



(6.8) 
(6.9) 



2 For a definition of spin- weight, see Appendix [h| 
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Spin- weighted spherical harmonics are defined in terms of the matrix representations of the rotation group, D l _ sm , 

by 



s Y lm (6, ¥>) = 7(2* + l)/^D l _ sm {4>, 9, 0), (6.10) 

where s is the spin weight and oYi m are the ordinary spherical harmonics Yj m . See Ref. j6f| for a definition and 
a detailed list of the properties of the quantities D_ am . For each value of the spin weight s, these spin- weighted 
spherical harmonics form a complete orthonormal set, i.e. 



iYi m * a Yi' m i dtt = 5i>i5 rn > m . (6-11) 
Thus we can expand ??o 5 ^i 5 ^-i on the spin-weighted spherical harmonics bases as 

MrAv) =E/o A M oY A (0,<p) 

A 

tfiM.vO =£/+(»•) iWp) 

A 

0-iM>¥>) =X) /-W -lW.p), (6-12) 

A 

where the subscripts + and — are understood to mean +1 and —1 respectively, and A denotes (Im). Therefore the 
expansion of the mode function can be written as 

£ = ^[/o m (r)olW +f l r(r) l Y lm ih + f lm (r) ^Y lm m] . (6.13) 
I'm 

In Appendix || we relate the coefficients f l s m {r) appearin g in this e xpans ion to the coefficients Ui m (r), Vi m (r), and 



Wi m (r) appearing in the expansion (p.l|) above [see Eqs. (Hll) and (H12)] 



We will be using a nonstandard notational convention involving the symbol A which we now explain. In rotating 
stars, I is generally not a good "quantum number" while m is. In general there will be other good quantum numbers 
that will allow us to enumerate the modes. What is meant by the symbol A is a set of numbers representing both 
the good and bad quantum numbers of the object in question. For example, when A is used on a s Yi m we mean 
the string of quantum numbers (l,m). However, when A is used on the mode function / A it may mean a string of 
numbers that includes, but is not limited to I and m. When we sum over A it is implied that we sum only over bad 
quantum numbers. Thus we might label the above component $o with another index (A say) representing only the 
good quantum numbers that are left after the bad ones have been summed over, e.g. 

#£(r,e,<p) = ]T/ A (r) o^aCM- ( 6 - 14 ) 



And finally, when we use the symbol A not in a sub- or super- script it is defined to be A = + 1). 

We now express the covariant derivatives of the mode functions that appear in the coupling coefficient (4.20) in 
terms of our new expansion coefficients f^(r). The details of the calculation are given in Appendix [l|. We define the 
functionals G A [/ A (r)], F s A [/ s A (r)], and i/ A [/ s A (r)] to be 



ff A [/ A (r)] = --±-/ A ( r) yA^2. (6.15) 



Then the components in the basis (1, m, m) of the tensor £^ can be written as 
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K ;l 


A 


S ;/ 


A 


pm 
> ;m 


= E G +[/] O^A, 

A 
A 


£171 

^ ;m 


= E G -[/] O^A, 

A 


£771 

;m 


A 


cm 
? ;m 


= E^-[/] -2^, 
A 


S ;m 




- 

S ;m 


= E^t/] -^a, 



(6.16) 

where we abbreviate G A [/ A (r)] as G A [/], etc. Note that this formalism allows us to write down forms for the covariant 
derivatives in terms of functions of r only. The angular derivatives have all been computed. Thus the formalism avoids 
the large curvature terms that can arise from angular covariant derivatives in spherical coordinates. Another advantage 
is that the formalism allows us to integrate the coupling coefficients over angles immediately; no integrations by parts 
are needed. To see how this works we note that the integral over angles of three of the spin-weighted spherical 
harmonics can be written simply in terms of Wigner 3 — j symbols fl65| , which are easily calculated: 

si Yi imi S2 Yi 2 m 2 s 3 ^3Tri3 dn = ( Sl Yi imi S2 Yi 27n2 S3 Yi 3m3 ) 



Ait^(2h + 1)/4tV(2^2 + 1)/4V(2^3 + 1)/4tt 

x / D l i Simi (^e : o)D l i S2m2 ( Vl e,o)D l i S3m3 (^e,o)dn 



= 4tt v /(2; 1 + l)/4vr v /(2; 2 + l)/47iV (2Z 3 + 1)/4tt 
h h h \ l h h h 



-81 -8 2 -S3 J \ m l m 2 m 3 

Note that this integral vanishes unless m\ + m 2 + m 3 = 0, cf. the m selection rule ( |5.12 ) 



(6.17) 



C. Integrating the coupling coefficients over the sphere 

To begin integrating the expression ( 4.2CQ over angles we write 



kabc = j r dr KABc(r), (6.18) 

and note that in calculating kabc we can ignore mode indices as long as we symmetrize at the end of the computation. 
To this end define an operator S that symmetrizes over mode indices: 

SL\ 1 a 2 a 3 = ^- , (A 1 A 2 A 3 ) = g (£AiA 2 A 3 + ^A 2 AiA 3 



Now expand the scalar V • £ A as 



and the scalar S<f>A as 



+ L\ 3 \ 2 A 1 + LA1A3A2 + £a 3 AiA 2 + iA 2 A 3 Ai) ■ (6.19) 
A 

a = 50 A = J2 a o( r ) oY A . (6.21) 

A 



The second covariant derivatives of this scalar are then given by 



a- 



inn I 



E ^o A [a] 2*A, a irnrll = E oi a ] -2>A 



A A 

a-i m = a.„ a = ^C A [a] iY\, a ;te = a ;f nl = E^-^ ~ iYa 
a A 



a-mfn = a-rnm = ^ Dq [a] Y A , CL-U = ^ a A rr Q Y\, (6.22) 
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where the functionals E^[a(r)], C A [a(r)] and Z) A [a(r)] are 

A 



£„ A K r )]E-VAM fl :( r ), 



A A 2 



sA (\ 



r \ t 



A _ A 



(6.23) 



We are now in a position to write down kabc{t)'- 

A1A2A3 

+65 [L^ AaAa (r)( y Al 1^-1^3) 

+65 [ii 3 ^ a 3 (r) ( y Al 2 fa 2 - 2 y A3 ) 

+65 

+ (+ 



^AiWOUyvayv.iyo 



-)}. 



(6.24) 



where the (+ <-> — ) symbol means repeat all terms in the above equation replacing s — > — s in the integrals and 
interchanging + and — subscripts in the functions (detailed below). The functions i A j AaA (r) are given by 



ainpo/ 50 5 ° 9 ° 2 P dr^° /o /o 



5 #1-1)3 



1 



Ai f A 



/o V2 .,/o A3 ., + 5o Al G^[/]G^[/] 



+2p 



1 



/o Al X 2 X 3 ,r + G+ [/]G+ 2 [/]Gt 3 [/] 



Ai r fi ^iA 2 



4?A 3 A 3 W=f(r 1 -iK i n 2 [/]/- 3 , r 
2 

cAi 1 ^<Air-fi 77iA 2 r f 1 j-A. 



Ai P A 2 



(6.25) 



-p ( ^ Al H/+ 2 /- 3 + fo'C+W- 1 



-p 



/o Al , r + G + 1 [/]J^+ [/]/- 3 , r 

Ai f A 2 fA 3 



L*gA 2 A 3 (r) = >(ri - l) 5o Al ^ 2 [/]ff- a [/] +pG + M/]ff+ 2 [/]J?- 3 [/], 



,Ai W A 



-(+ 



^A 3 A 3 W =PH^[f]tt\ r K 3 [f}-P~K 1 W+ 2 tf 3 



(6.26) 

(6.27) 
(6.28) 



where a comma indicates ordinary differentiation. 

To summarize, given a set of three modes £ Al £ c , the computational method is to (i) Use Eqs. (Hll) and ( H12| ) 
to compute the expansion coefficients /j m (r) for each mode in terms of the more standard expansion coefficients 
Ui m (r ), V/ m (r), and Wi m {r) appearing in the expansion fl6.1| ), and (ii) Use the expressions ( |6.17 ), (6.18) and (6.24) - 
( |6.28 ) to compute the coupling coefficient kabc from the functions f l s m (r). 



VII. SUMMARY 

In this paper, we have formulated a perturbative approach to the nonlinear interactions of unstable r— modes in 
a neutron star. Our formalism presumes that mode growth saturates at moderately low amplitudes, so that we can 
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model the modal interactions via three-mode couplings. By developing further a previous perturbation theory for 
rotating stars, we have found equations of motion for the mode amplitudes of rotating stars that are uncoupled at 
linear order when acted upon by an external force. This feature of our formalism is essential for following the cascade 
of energy from one mode to another when lowest order nonlinear couplings are included. 

The important astrophysical question is what determines the saturation amplitude for unstable r— modes. Fun- 
damentally, if mode-mode coupling is the dominant damping, we expect the saturation amplitude to be set by a 
competition between the growth rate of the instability, and the amplitude-dependent rate of drainage of energy from 
an unstable mode to other stellar modes. Our formalism gives explicit formulae for the lowest-order coupling coef- 



ficients among stellar modes. The numerical simulations by Stergioulas and Font |16|, and Lindblom, Tohline and 
Vallisneri [[l7|, have found the coupling of the unstable I = \m\ = 2 r— mode to other stellar modes to be surprisingly 
weak, resulting in a timescale for energy transfer to other modes that is at least ~ 20 rotation periods even when the 
modal amplitude is substantial. 

Although we shall present detailed numerical solutions of our equations of motion for modal amplitudes in a 
subsequent publication, there are already hints, in the formal developments presented here, of an explanation for the 
apparent weakness of the coupling of the I = \m\ = 2 i — mode to other stellar modes. Strong coupling requires either 
near-resonance or a large coupling coefficient, or both. Although the rotational modes have frequencies comparable to 
the r-mode, we have found that the possible couplings to such modes may be limited. For example, in zero-buoyancy 
stars parity arguments prevent the coupling of three r— modes, and the couplings involving two r modes and one 
hybrid rotational mode vanish to zeroth order in the stellar angular velocity. In nonzero- buoyancy stars, the coupling 
between r-modes is again vanishing to zeroth order in f2, but is nonzero when correction terms of order 0(fl 2 ) are 
included |64|. Thus, the couplings among r— modes are small, generically. The coupling of two r— modes to a /—mode 
need not be small. However, since the eigenfrequencies of /—modes are ~ (GM/i? 3 ) 1 / 2 , which is much larger than 
for a slowly rotating star, the excitation of an /—mode from small amplitude by an unstable r— mode is likely to be 
suppressed, by factors ~ Q 2 R 3 /GM. 

Finally, as explained in the introduction, the possibility of the r-modes saturating at small amplitude is not neces- 
sarily incompatible with the numerical simulations to date, or with the apparent weakness of the nonlinear couplings 
of the r-modes. If the modes do saturate at small or moderate amplitude, then the formalism developed in this paper 
may be sufficient to explore the saturation process. On the other hand, it is possible that the formalism developed 
here may have to be developed further in order to understand the saturation. Some possibilities are 

• It may be that the r-modes stop growing only in the very nonlinear regime, where shocks develop, as in the 
simulations of Lindblom, Tohline and Vallisneri |l7| . If this is the case then the perturbative formalism developed 
here will not be useful. However, as explained in the introduction, the simulations to date do not show that the 
strongly nonlinear regime is reached. 

• It may be that an analysis with only three- mode couplings will be insufficient, as cubic potentials generically 
have instabilities at large amplitude pq| , but that including four- mode couplings will be sufficient to allow an 
exploration of the saturation process. 

• This paper has only studied the coupling coefficients to zeroth order in the stellar angular velocity. However, it 
may be that the dominant energy-transfer channels involve coupling coefficients that are nonzero only at 0(fl 2 ). 

• We have focused attention in this paper on zero-buoyancy stars, for simplicity. However the quantitative details 
of the saturation process in real neutron stars may be altered by the presence of buoyancy forces. 

• Finally, we have used Newtonian gravity throughout. It is conceivable that important energy transfer channels 
could come about via coupling coefficients which vanish at Newtonian order (to zeroth order in f2), but which 
are non-vanishing when post-l-Newtonian corrections are included. This will depend on the relative sizes of the 
two dimensionless parameters GM/(c 2 R) and R 3 fl 2 / (GM). The fact that the gravitational radiation reaction 
force on the I = m = 2 r-mode is dominated not by the "Newtonian" quadrupole coupling but instead by the 
"post-Newtonian" gravitomagnetic coupling is a hint in this direction. 
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APPENDIX A: HAMILTONIAN ANALYSIS OF LINEARIZED PERTURBATIONS 



The purpose of this appendix is to derive the mode decomposition formalism described in Sec. [n] of the body of the 
pap er. As discussed in Sec. 0, most of that formalism is contained in a series of papers by Sc hutz and collaborators 
p3| ^6| . Schutz introduces the phase space mode expansion, shows that the set of vectors ( |2.3l| ) form a basis for 
Ti x H, and gives an extensive discussion of Jorda n ch ains. The main new feature that we introduce is the explicit 
computation of left eigenvectors [defined in Eq. (A.31) below] in terms of right eigenvectors, and the subsequent 
derivation of the equations of motion in the explicit form (2.37) abo v e. Fo r completeness, we sketch in this Appendix 



derivations of all the building blocks of the formalism. Sees. Al - A 4 detail the formalism for non-Jordan chain 



modes. However, evolving the differential rotation and/or total spin of the star requires an ex tens ion of the formalism 
to compute the equations of motion for Jordan chain modes. That extension is given in Sec. A 5. 



1. Phase space equations of motion 



We start by reformulating the equation of motion (2.18) as a pair of first order equations rather than a single second 
order equation. Our treatment closely follows |24[] , except that we use canonically conjugate variables instead of £ 
and £. 



The equation of motion (2.18) can be derived from the Lagrangian density 



C ■ - 2 i-\i-i - 2 i-C-C ■ a ext (t) • £• 



The momentum canonically conjugate to £ is 



8{ 2 s ' 



and the associated Hamiltonian density is 

f 



§1* 2^ 



1 



C C 4 a, s: ;/:•£. 



The Hamiltonian equations of motion can be written as 

C = T-C + F(t), 

where 



the operator T is 



and where 



C(*,x) 



-C 



F(t) = 



7r(t,x) 



B 

-IB 2 



1 





a ex t(i) 



If we now specialize to the case of no forcing term, a ext (i) = 0, and assume a solution of the form 



C(*,x) 



f C(x) 



then we get from Eq. ( A4) the eigenvalue equation 

[T + iu>] ■ C(x) = 0. 



(Al) 



(A2) 



(A3) 



(A4) 



(A5) 



(A6) 



(A7) 



(A8) 



(A9) 



It is straig htforw ard to show using Eqs. ( |A2| ) and (|A5|) that Eq. (A4) is equivalent to the configuration space equa tion 
of motion (2.18), and that the eigenvalue equation (|A9|) is equivalent to the quadratic eigenvalue equation (2.18). 
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2. Right and left eigenvectors and Jordan chains 



Label the distinct right eigenvectors of T as Cai an d the associated eigenfrequencies as uja, so that 

[T + iuja] ■ Ca = 0. (A10) 

Since the operator T is not Hermitian, we will have also left eigenvectors Xa (distinct from the right eigenvectors) 
that satisfy 

[Tt - icj* A ] ■ xa = 0. (All) 

Here 



±B -C + iB 2 
1 5 B 



(A12) 



is the Hermitian conjugate of T, since C T = C and B^ = -B. 

Since the operator T is not Hermitian, the set of its right eigenvectors will not in general be a complete basis. 
However, one can obtain a complete basis if one includes all Jordan chains |2q] . What this means is as follows. For a 
given eigenvalue — iu>A, let Va be the subspace of Ti.2 = Ti. © Ti. consisting of vectors C that satisfy 

[T + iu A ] m ■ C,=Q (A13) 

for some integer m > 1. Clearly the space Va contains all the right eigenvectors associated with u A . Now it can be 
shown that the direct sum of all the subspaces Va (one for each distinct eigenvalue) gives the entire space Ti.2- Hence, 
to obtain a basis for Ti.2, it suffices to find a basis for each space Va- One can always find such a basis of the form 
{Ca,o-} where < a < pa, Ca o = Ca I s the right eigenvector, Ca,i> • • • > Ca,va are a se ^ vec tors in Va that satisfy 

[T + iUA]-CA,* = <A,*-l, (A14) 

for 1 < a < pa , and that form the Jordan chain of length pa associated with £ A . For each eigenvector £ A , either pa = 
and there is no associated Jordan chain, or there is a chain of length pa > 1. The set {Ca.o- A = 1, 2, 3, ... , < 
v < Pa} of Jordan chains, including the right eigenvectors, forms a basis of 7^2- 

These Jordan chains are right Jordan chains. We also have for each A a left eigenvector Xa = Xa o as discussed 
above, and left Jordan chains consisting of vectors Xa a f° r 1 — a — PA satisfying 

[T f - iu* A ] ■ Xa,* = Xa,*-i, (A15) 

for 1 < A < pa- Note that the length pa of the left Jordan chain must be the same as that of the right Jordan chain, 
for each A. The basis {xa a-} °f left Jordan chains can be chosen to be dual to the basis of right Jordan chains in the 
sense that 

{XAm i Cb,\) = S ABS pA ,a+\; (A16) 

for all A, B and for < a < pa and < A < ps Pq] . The inner product here is defined in the obvious way as 

<K,t] , [£',7r']> = {£,0 + <*->*■'>, (A17) 
where the inner products on the right hand side are given by the definition ( ^20| ). 

3. General mode expansion and equations of motion for mode coefficients 

We expand C(^: x ) on the basis of right Jordan chains as 

PA 

C(i,x)=^^c A)(T («)C A , (T (x), (A18) 

A cr = 



where from the orthogonality relation (A16) we have 



29 



C A,er — (XA,PA-<T ' C) ' 



(A19) 



Subst ituting the expansion ( A18 ) into the dynamical equation (|A4|), and using the defining properties ( A10 ) and 
(A14) of right Jordan chains yields 



(A20) 



A <r=0 



A a=0 



where we define £ A _ 1 = 0. Now multiplying on the left by Xb a, using the orthogonality relation (A16), and then 
relabeling the indices yields the equations 



for < a < pa — 1, and 



C A ,a + iUACA,a - C A ,a+l = (XA,p A -a > F (*)) > 



c A , PA + iuj A c A , PA = (xa,o > p (*)) 



(A21) 



(A22) 



Thus, the equations for the mode coefficients CAfi{t), CA,i(t), . ■ ■, ca .va (t) are c oupled together for each fixed A, but 
the different A's are uncoupled. The general solution of Eqs. ( |A2l| ) and ( |A22| ) in the case of no forcing terms is a 
polynomial in time multiplied by the usual complex exponential: 



VA-a 



CA, 



(A23) 



A=0 



where 7,4,0, Ja.i, ■ ■ ■ , lA,p A are constants of integration. 



4. Specialization to no Jordan chains 



We now specialize to the situation where there are no Jordan chains. In this case the orthonormality relation (ALE) 
reduces to 



and the mode expansion ( A18) becomes 



(Xa > Cb) = &AB j 



C(t,x)=^c A (t)C A (x). 



(A24) 



(A25) 



The inverse of this mode expansion is 

CA(t) = {XA,<(t)), 

and the equation of motions (A21) and ( A22| ) for the mode coefficients reduce to 

c A + iuJACA = (Xa i F W) ■ 



(A26) 
(A27) 



a. Translation of results to configuration space variables 



Our goal now is to write these results for no Jordan chains entirely in terms of the Lagrangian displacement £ and 
its time derivative £, and the modes (£ Ai lua) of the configuration space formalism. To do this we proceed as follows. 
We write the right eigenvector £ A as 



Za 

™A 



(A28) 



and similarly write the left eigenvector Xa as 
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Xa 



cta 

T A 



(A29) 



Then the definitions (A10) -( [All ) of right and left eigenvectors together with the formula (A6) for the o perat or T 
s hows that £ A and ta are right and left eigenvectors of the operator L(ua) — —w\ — zw^B + C [cf. Eqs. ( 2.28 ) and 
p~TTTi| above] 



L(wa) • £ A 

L(wA) f ■ TA 




0. 



In addition we obtain the formulae 



and 



it A = -iu) A $A + B ' €a/ 2 



a- a = iu*A T A - B ■ Ta/2, 



(A30) 
(A31) 

(A32) 

(A33) 



which using Eqs. (A34) and (A35) allow us to write the phase space right and left eigenvectors C, A an d Xa entirely in 
terms of the configuration space left and right eigenvectors £ A and t a : 



and 



Ca = 



Xa 



£a 



iu* A T A - |B • TA 
TA 



(A34) 



(A35) 



Now using the definitions (AS) and flAq) and the formula (A34) we can rewrite the mode expansion (|A25|) as 



kit) + ±B • £(t) 



-iu A £ A + h B ■ Ca 



(A36) 



which is equivalent to the expansion (2.34) quoted in the body of the paper. Note that there is a one-to-one corre- 
spondence between right eigenvectors solutions £ A of Eq. ( Ald|) and s olutio ns (£ A ,u> A ) of Eq. ( 2.18 ), which allows us 



to identify the sums over A that ap pear in the expa nsions (2.34) and (A36) 



Next, w e can use the formulae (A34) and ( A35 ) to rewrite the orthogonality relation ( A24 ), the inverse mode 
expansion (A26) and the equation of motion (A27). The results are 



Sab = (ta , B ■ £ B ) - i(uj A + ujb) (ta , £ B ) 



and 



c A (t) = -i (ta , w A Z{t) + *B • C(t) + 
c A + iuj A c A = (t a , aext(i)) , 



(A37) 
(A38) 

(A39) 



where we have also used Eqs. (2.22), (A2), (A5) and (A7) 



b. Computation of the left eigenmodes 

The last step in the construction is the explicit computation of the left eigenvectors t a in terms of the right 
eigenvectors £ A . For a giv en eigenfrequency u, let Vr(uj) denote the space of associated right eigenvectors £ satisfying 
L(w) • £ = [cf. Eq. (|t2f)], and let Vl(w) denote the space of associated left eigenvectors r satisfying L(w)^ • r = 0. 
The spaces Vr(lo a ) and Vl(^>a) have the same dimension, which is the degeneracy associated with the eigenfrequency 
lj a . However, Vr(lu a ) and Vl(^a) will not in general coincide for complex oja- We now show that for real uja, the 
spaces Vl(wa) and Vr(oja) do coincide, which means that the left eigenvectors associated with a given eigenfrequency 
can be expressed as linear combinations of the corresponding right eigenvectors. This is the key result that we use to 
compute the left eigenvectors. 

We shall need the following properties of the spaces Vr(ui) and Vl(w) of right and left eigenvectors discussed by 
Schutz |^5| . Suppose that £ lies in Vr[uj). Then 
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The function £ also lies in Vl(u* ), ie, is a left eigenvector with eigenfrequency w*. This can be derived by taking 
the Hermitian conjugate of Eq. (|2.18 ), and by using the fact that B is anti-Hermitian and C is Hermitian, which 
yields L(w*)t ■ | = 0. 



The function £(x)* lies in Vr(— ui*). This follows directly from the fact that the equation of motion (2.7) is real 



• Let g be the mapping that takes (r, 9, ip) to (r, 0, —ip) in spherical polar coordinates, or equivalently (x, y, z) — > 
(x, —y, z) in Cartesian coordinates, i.e., reflection in the xz plane. We define g*£ be the pullback of £ under g 



so that if 



£ = f (r, 6,<p) A + e (r,e,<p) ^ + ^ (r, *, ¥>) ^ , 



(A40) 



then 



<96» 



(A41) 



Since the background star is invariant under the transformation tp — ► — tp, t 

v a (- W ) pa. 



— i, it follows that lies in 



• By combining the previous properties, it follows that (g*£)* belongs to both Vr(w*) and Vl(w). The mapping 
£ — » (<7*£)* was nrs t written down by Schutz in Ref. |23] where it was denoted S 1 . 



We next recall the notation used in Sec. II C 3| above. We write the distinct eigenfrequencies as u> a , and the right 
eigenvectors as £ A = £ a k where 1 < k < n a and n a is the degeneracy associated with the eigenfrequency u) a . 

Consider first the case of complex uj a . In this case it follows from the above that the mapping £ — > (<7*£)* maps 
Vii(wo) onto Vi(ijJa) PH]- Hence we can write 



T a ,k 



k' = l 



Tkk> {g*€a,. 



(A42) 



for some n a x n a matrix T^k' ■ By substituting Eq. ( |A42| ) into the orthogonality relation (A37) we can solve for the 
matrix T^k 1 and thus obtain the left eigenvectors T a ^- We note that Schutz uses a differe nt ba sis of left eigenmodes 
for which T a ,k — (g*£ a k)* [Eq- (3-14) of p5fl ], but for which the orthonormality condition (A37) do es n ot hold. This 



method of obtaining the left eigenmodes for complex frequencies is discussed in more detail in Sec. A 5 



below. 



For real frequencies, the situation is much simpler. Every £ in Vr{oj) also lies in Vl(w*) = Vl(u>), and hence the 
spaces Vr(ui) and Vl(w) coincide 0. Hence we can write 



a,k = 



L0 kk' Sa 



k' ■ 



(A43) 



for some n a x n a mat rix fiS- (The factor of — i and the complex conjugation are inclu ded for convenience.) By 



substituting Eq. (A43) with t a 



T a,k 



and £ B = £ fe ; into the orthogonality relation ( A37) we obtain 



Sab Ski = B { kk> M a 



ik' ,bl i 



(A44) 



where the matrix M. a k' ,bi is defined in Eq. ( 2.4C| ). It follows from Eq. ( A44 ) that, first, the matrix Ai a k,u is block 



diagonal in the sense claimed in Sec. II C 3, and second, for a given value of a the matrix Bffi is just the inverse of 
the the diagonal block M a k ai- 

Now as explained in Sec. II C 3, we can always choose the basis £ a k of Vfl(w a ) to diagonalize the matrix M. a w ,ai- 



In particular, this will be automatically true when n a = 1 and there is no degeneracy. For such bases, the matrix 
Bkk' will also be diagonal. It then follows from Eq. ( A43) that each left eigenvector can be written as 



23 For a mode (£ A ,u> A ) with no degeneracy, it follows from this that we can choose the phase of £ A so that (<?*£ 4 )* = £a- 
Hence, if m is the azimuthal quantum number and we write £ A = exp[imip]£ A , the r and components of £ A are purely real 
and the tp component is purely imaginary. 
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TA 



b A 



(A45) 



for some constant b A . Su bstituting this into the orthonormality relation (A37) yields the form ula ( 2.36) for b A . 
Finally, the form ula (|A45 ) for left eigenvectors can be combined with the ortho g onali ty rela tion ( A37 ), the inverse 



mode expansion ( A3q ) and equation of motion (A39) to yield the versions (2.32), (2.35) and (2.37) of these relations 
quoted in the body of the p aper Pj . 

Finally, we note that Eq. (A44) shows that the matrix A4 a k .ai is non-degenerate. Thus, if that matrix is degenerate, 
then the assumption underlying the derivation of Eq. ( A44 ) — that there is no Jordan chain associated wit h the 
eigenfrequency lu (1 — must fail. In particular, for non-degenerate modes, if the constant b A defined by Eq. ( 2.3(: ) 
vanishes, then the mode must be a Jordan chain mode. 



5. Computation of the left eigenmodes for Jordan chain modes 



In this subsection we show how to generalize the above analysis to obtain the left eigenvectors x A a f° r Jordan 
chain modes. The construction is useful, for example, in computing the evolution of stellar differential rotation, which 
is described by an infinite set of Jordan chains of length 1 (see Appendix ^) . 

The essential idea is to use the phase space variables rather than the more complicated configuration space variables, 
and to use an operator which maps right eigenvectors onto left eigenvectors. We define the operator 



M 



1 

-1 



(A46) 



which is essentially the symplectic structure defined by Friedman and Schutz It satisfies ]VF = M = M 1 and 

• T M = -T f . (A47) 

It follows from Eqs. ( A47 ), ( A10|) and (All) that if £ is a right eigenvector of T with eigenfrequency u>, then M • £ is 
a left eigenvector with eigenfrequency lj* . 

Now the set of vectors M • C, A a forms a basis, since M is invertible and the set of vectors C, A a is & basis. Therefore 
we can write the left eigenvectors and left Jordan chain vectors as 



X A , a =J2 B *A^BxM-C 



(A48) 



B,X 



for some matrix B A(7 ,b t, w here the complex conjugation is included for later convenience. Inserting this into the 
orthogonality relation (Al£) gives 



where 



B A <t,Ct M-Ct,BX — S A B &a+X,p A , 
Mct.BX = (Cc,t,M-C B , X ) 



(A49) 



(A50) 



are the matrix elements of the operator M. Equation (A4E) says that the matrices B and M are inverses of each 
other (up to the index per muta tion er — > p A — er). 

The reason the ansatz (A48) is useful is that the matrix M. is almost diagonal, in the sense that it satisfies the 
identity 



(u>* A - uj B )M ArT . 



BX 



0. 



(A51) 



It follows from Eq. (A51) that M. is block diagonal, with one block for each real frequency and one block for each 
pair (w,w*) of complex frequencies. Therefore, we can obtain the left eigenvectors by inverting each diagonal block 
of M Aa ,BX to obtain B Aa ,BX, and by using Eq. (A48). 



24 The validity of th e ort hogonality relation (2.32), for those values of A and B for which uja <+>b, can also be derived directly 
by contracting Eq. (2.30) on the left with £ B , subtracting from this the same equation complex conjugated with A and B 
interchanged, and dividing by uja — ub- 
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To derive the identity (A51) we use Eqs. (A47) and (A14) to show that the vectors M ■ Ca a form a left Jordan 
chain with frequency lo* a : 



[Tt - iu A ] ■ (M • C A ,a) = — M • C 



A.cr-l- 



(A52) 



Next, we have 



J A M Au , B \ = u A (C A ^ , M-C B ,x) 

= (u JA M-C A ,a , C B ,x) 



(-iM-C 



A,a-1 



/'I' • M • <' , „ , Cb,x) 



= i {C A ,a-i , M • Cb.x) + < (M • Ca,* , T ■ Cb,a) 
= * (Ca,*-i , M • Cb.a) + * (M • C A , CT , -*wbCb,a + 
where we have used Eqs. (|A14[ ) and ( |A5^ ). The result ( |A5^ ) can be written as 

(W^ - UJ B )M A a,BX = i (M A (a-l),B\ + M Aa>B (X-l)) , 



,B,X 



-l) 



(A53) 
(A54) 









TT 







and the identity (A51) follows from iterating the identity (A54). 

Schutz suggested a different general method of constructing left Jordan chains. That method is based on the 
anti-linear operator S2 defined by 

(A55) 

which maps right Jordan chains Ca a of T with frequency uja onto right Jordan chains S2 • Ca a with frequency u> A . 
Therefore composing this map with the operator M yields a mapping C — > M • S2 • C which takes right Jordan chains 
of frequency u> to left Jordan chains of frequency lo. The approach here omits the mapping S2 and is simpler to use 
when the eigenfrequencies are real. 

For Jordan chains of length zero (ordinary eigenvectors), the method (A48) reduces to the method (A43) of Sec. 
A 4 above. Dropping the Jordan ch ain in d ices a and A, th e ansatz ( A48| ) reduces to Xa = ^2b^ab^ ' Cs with 
BacMcb = Sab, which using Eqs. dA34[ ), ( |A35D and flA46| ) implies 



T A = -i^B* AB i B , 



cf. Eq. (A43) above. Similarly, the definition ( |A5C| ) can be simplified using Eqs. (A34) and (A46) to give 

Mab = (CaM-Cb) 

= (£ A , zB • £ B ) + {cj* a + uj b ) {£ a , £ B ) , 



(A56) 



(A57) 



which is a generalization of the definition (2.40) valid for complex frequencies. 

As an example, consider now the case of a mode (Ci> w i) with complex eigenfrequency which is non-degenerate. 
Then there is an associated mode {^2^2) with LO2 = w*, and by choosing the normalization of £ 2 we can without loss 
of generality take 

£2 = (9*Si)*, (A58) 
cf. Sec. A4b above. The corresponding 2x2 block of the matrix M. is from Eqs. ( |A51 ) and (A57) of the form 



M 



AB 



V 

V* 



where V = M.\i = N[^ 1 ,u>i], where the functional N is 

N[£, w] = <£ , »B ■ (<?*£)*) + 2cj* (£ , (g m & 
The corresponding left eigenmodes are therefore 



and 



Tl -~v i2 



T 2 



I 

V 



■ix- 



.£1) 



(A59) 
(A60) 
(A61) 
(A62) 
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a. Zero frequency Jordan chains of length one 



We now carry through the explicit computation of the left Jordan chains for the case that arises in practice in stable 
rotating stars, that of zero- frequ ency Jordan chains of length one. Consider first the case of no degeneracy. Dropping 
the index A, we have from Eq. (A14) two Jordan chain vectors Co an d Ci which satisfy 



T-C o = 0, T-d=Co- 

The 2x2 matrix M aT — (C CT > M • £ T ) satisfies the identity 



(r-1) 



Q 



from Eq. ( A54 ). Using this identity together with the fact that M is Hermitian gives 

Mar = 



i/3 
-i/3 7 



(A63) 
(A64) 

(A65) 



where (3 and 7 are real. Under the transformation Ci — > Ci + xCoi which preserves the defining relations (A62), we 
have 



P~>p, 7^ 7 + 2/3Im(x). 



(A66) 



It follows that we can choose Ci to make 7 = 0. Now combining Eqs. jA48|) , ( A49 ) and ( A65 ) gives for the left Jordan 
chain vectors 



Xo = -p M • Co 
Xx = ~M. Cl . 

We now write these results in terms of configuration space variables, using the notation 



7T n 



(A67) 
(A68) 

(A69) 



0"n 



for a = 0, 1. First, the relations (A63) can be written as 

C-£ = 0, C.€i = -B-€ 0J Wo=B.f /2, 7n=€o + B-^/2. 



Second, the formula for (3 is, from Eqs. (A46), (A69) and (A71) 



(3 = -1 (C , M • d) 
= (£0 , €o> + «o . B ■ €1) 



Third, the relations (A67) and (A68) together with Eqs. (A46), (A69) and (A70) imply that 



1 c 1 t 



The corresponding equations of motion are, from Eqs. (A21) and (A22), 



Ci — -p (Co > a eirt) 

CO = ~o (Cl , Bext) + C\. 



(A70) 
(A71) 

(A72) 
(A73) 

(A74) 
(A75) 
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We now generalize these results to allow for degeneracy, since the space of zero-frequency Jordan chain modes in 
rotating stars is highly degenerate (see Appendix |j) . We use the notation 



A = (a, k) 



(A76) 



of Sec. II C S above, where a labels the distinct eigenfrequencies and k the eigenvectors associated with each eigenfre- 
quency. Therefore the Jordan chain vectors can be written as 



C Aa Caka 



(A77) 



where a labels the Jordan chain vectors associated with each eigenvector. For simplicity we drop the label a in what 
follows. The defining relations for Jordan chains of length one are 



T-C fco =0, T-C kl =Cw 
The matrix Mka,ir — (Cfe CT : M • £ lr ) satisfies the identity 



■A4fc( CT -l)jT + MkaMr-l 



(r-1) 







from Eq. (A54). Using this identity together with the fact that M. is Hermitian gives 



Mko,n 

Mki M 



= 



where the matrices (3ki and jkl are Hermitian. 
Next, the transformation 



CfeO 



-Mk\,lo = iflkl 
Ikh 



CfcO — FklClO 

Cki = -ffczCii 



(A78) 
(A79) 



(A80) 
(A81) 
(A82) 



(A83) 



preserves the relations ( |A78| ) . Under this transformation the matrix (3^1 transforms as f3 — > f3' = P • (3 ■ F^, and it 
follows that we can choose the basis Cko to diagonalize (3m, so that 



Similarly the transformation 



Pkl — fik&kl- 



CkO ~ * CfcO — CfeO 

Cfei — * Cki = Cki + FkiCm 



preserves the relations (A78). Under this transformation the matrices transform as 

7 -> 7' = 7 + iF ■ (3 - i(3 ■ F f . 



(A84) 



(A85) 



(A86) 



Now it follows from Eqs. ( A80| )-( [A82 ) and the fact that the matrix Mka,w is invertible that (iu is invertible. Therefore 
we can find a transformation matrix F^i which achieves 7' = in Eq. (A86), namely 



F = i7-/3"V 2 - 



(A87) 



Therefore we can choose the basis Cki to make 7/^ = 0. The rest of the analysis now proceeds exactly as for the 
non-degenerate case above, and we obtain the equations of motion 



1 

Jk 



(CfcO > a exty 



fa 



(£fci j a cxt ) + Cfci. 



(A88) 
(A89) 
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APPENDIX B: LINEAR INDEPENDENCE OF A SUBSET OF THE RIGHT EIGENVECTORS 



In this appendix we show that the of non- Jordan-chain right eigenmo des £ A for which the constant 6 a is positive 
are linearly independent. We start by recalling the notation used in Sec. II C 4 above: the index a labels the distinct 
right eigenvectors £ A for which 6a > 0, and we write A — (a,e) where e takes on the values e = + and e = — . Then 
£ Q = £«,+ and £ Q - = C a , while w a ,± = ±u a and b a .± = ±b a . 



By combining the orthogonality relation ( A37 ) with the formula ( A45 ) for left eigenvectors ta, and choosing 
A = (a, +) and B = (/3, +) we obtain the following form of the orthogonality relation: 



(Bl) 



Suppose now that the vectors £ a are not linearly independent. Then, one can find coefficients c a , not all zero, so the 
vector 



<1 = H c « £<* 



(B2) 



is zero. Now contract Eq. (Bl) with c* a cp and sum over a and (3. The result is 

M 2fo " = (q, iB • q) +2/ q, ^c a u a £, a 

a \ a. I 

-0, 



(B3) 



where the second line follows from q = 0. Since on the left hand side b a > for all a, this forces c a = for all a, 
which contradicts our assumption above. 



APPENDIX C: ALTERNATIVE FORM OF EQUATIONS OF MOTION 



In this appendix we present a version of the general equations of motion which is second order in time and which is 
more similar in form to the standard equations of motion for non-rotating stars. We start with the equations derived 
in Sec. II C 4 for the coefficients c a ,+ (t) and c Qj _(i): 



c Q ,+ + iw a c a . + = r-fa,+(t), 

Or* 



lUJryC. 



(Cl) 
(C2) 



where f a ,±(t) = (£ a ± , £ie X t(t)). As noted in Sec. II C 4, the Lagrangian displacement £(x, t) will be real if and only 
if c Qj + = c* _ and / Qj + = /^_. Here, however, we allow arbitrary complex £(x, t). We define 



e a (t) = c a . + (t) + c a -(t) 

Oa Oa 



(C3) 
(C4) 



where we have used Eqs . (|2.45| ) and (2.46). It is now straightforward to show using the definition ( |C3j ) that the 
equations of motion (|Cl|)-(|C2D are equivalent to 



e a (t) +u 2 a e a (t) = T a {t), 
which is the standard forced harmonic oscillator equation of motion. The forcing term here is 



3~ a 



% % 
■J- (fa,- + lUafa,-) + T~(fa,+ - i^afa, + ) 



2uj a 



(Re(| Q ) , a 



b a 



(Im(£J , a 



(C5) 

(C6) 
(C7) 



It is clear that e a (t) and !F a (t) will be real if and only if a ext (i) and ^(x, t) are real. We can express the Lagrangian 
displa cement £(x, t) in terms of e a (t) and e a (t) at any time t by using the first component of the mode expansion 
(2.44), and by using the relations 
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2 V cj„ 



1 



2 l6n 



'2uj a b ct 

i 

*2uj a b a 



(fa,+ ~ fa,-) , 
(/<*,+ — /a,-) • 



(C8) 
(C9) 



In this formalism, the number of "modes" e a (t) is the same for rotating stars as for non-rotating stars. In particular, 
for non-ro tating stars the expansion coefficients e a coincide with the standard expansion coefficients q a given by 
Eq. (2.25), when the mode functions £„(x) are chosen to be real. 

Note that the equation of motion (Cq) for a given mode has a characteristic feature which is peculiar to rotating 
stars, namely the forcing term depends on the time derivative a oxt (x, t) of the externally applied force per unit mass, 
as well as on a oxt (x, t) itself. For no n-ro tating stars one can choose a real mode basis £ a which removes the dependence 
on a ext [cf. the second term in Eq. (C7)], but for rotating stars it is not possible to find real mode bases. 



APPENDIX D: A SET OF ZERO-FREQUENCY, JORDAN-CHAIN MODES 



In this appendix we justify the claim made in Sec. II C 2 above that there are always Jordan-chain modes present in 
rotating zero-buoyancy stars, even for stable stars. We show that the well known class of purely axial, zero frequency 
modes that move the star to nearby equilibrium states with different angular velocity are Jordan-chain modes. Note 
that in this appendix we do not need to assume that the angular velocity of the star is small. For non-rotating 
stars, the existence of these Jordan chains was previously noted in a footnote by Ref. p£j 



It is most convenient to work directly with the linearized Eulerian equations of motion (2.5) and (2.6). Substitute 
into these equations the ansatz <5p(x, t) — 6p(x.) and 



<5u(x, t) = <5u(x) = r±5il(r±)e v 



(Dl) 



where (rj_, tp, z) are cylindrical coordinates with r± the distance from the rotation axis, and e v is the unit vector in 
the ip direction. One finds that this yields a time-independent solution, and that one can specify the perturbation 
5Q(r±) arbitrarily and solve for the density perturbation <5/?(x) 

Now switch back to the language of Lagrangian perturbation theory. The equations for a Jordan chain of length 
one with lua = are as follow s. [W e drop the index A and keep only the index a for convenience.] The phase space 
right eigenvector is, from Eq. (A34), 



Co 



and the associated right Jordan chain vector is 



Ci 



Here £ and ^ satisfy, from Eq. (A.14) 



and 



L(0) ■ £ = C • £ Q = 0, 



L(Q)-^ = C-^ = -B^ 



(D2) 



(D3) 



(D4) 



(D5) 



Now let £o( x ) be a mo de func tion of the form (Dl). Th en o ne finds that V • £ = and V • (p£o) = 0, an< ^ hence 
from the definition (2.13)— (2.15) of the operator C, Eq. (D4) is satisfied. One can then always solve Eq. ( p5| ) to 
determine £ 1; since the left hand side is a pure gradient for zero-buoyancy stars and the curl of the right hand side 



vanishes by Eq. (Dl). Thus the solutions (pl|) correspond to a Jordan chains of length on e. Th e vector ^ encodes 



the density perturbation. The general solution £(x, t) associated with £ is then, from Eq. ( A23 ) 



3 The density perturbation always vanishes if the background star is non-rotating but not in general. 
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£(x,i) = 7olo (x) +7i Ki(x) + €o(*)*] 



(D6) 



where 70 and 71 are constants of integration. 

The piece of the solution proportional to 70 has vanishing Eulerian density and velocity perturbations, and is pure 
gauge in the sense of Ref. [|2l|j . By contrast, the piece proportional to 71 is a physical, non-gauge mode. 

For non-rotating stars, the situation is a little differen t: a ll zero frequency modes are associated with Jordan chains 
of length one. This can be seen from the fact that Eq. (D5) always has solutions for non-rotating stars, since B = 0, 
whereas for rotating stars a solution only exists if the curl of B • £ vanishes. The associated solutions of Eq. (2.7) 
are of the form £(t ) = + t£ where £ and ^ satisfy C • £ = and C • ^ = 0. 

Finally, Ref. |26[ shows that in situations where there are Jordan-chain modes, if one adds any small perturbation 
to the system then unstable modes are generically created. Hence, Jordan-chain modes should generically not occur 
in stable stars. However, Schutz notes that a necessary condition for this argument to apply is that a certain matrix 



element of the perturbation not vanish. Specifically, if AT is the perturbation to the operator 
right Jordan chain with associated left Jordan chain \ A CT , then the argument requires 

(xa,o> AT-Ca,o)^0. 



V and Ca.v is a 



(D7) 



In the current context, the argument is evaded because of the fact that the matrix clement (D7) vanishes for arbitrary 
linear perturbations AB and AC to the operators B and C, from Eqs. (Ac) and (|D2]). Generic physical perturbations 
to the syst em do not correspond to generic mathematical perturbations to the operator T, because of the structure 
of Eq. (|A6|). 



APPENDIX E: THE ZERO FREQUENCY SUBSPACE FOR A NON-ROTATING STAR 

This appendix characterizes the space 7io of zero frequency modes for non- rotatin g s tars. The proofs in this 
appendix are Newtonian versions of the fully relativistic proofs given in Sec. Ill of Ref. |p7| . 
The Lagrangian displacements £(x) for zero frequency modes satisfy 



C£ = 0, 



(El) 



where the operator C is defined by Eq. (2.11) specialized to a non-rotating background star. From Eqs. (2.11)— (2.15) 
Eq. (El) can be written as 



V \gi(r)5p(x) + ty(x) + <? 2 (r)f»] - /3(x)g 2 (r)e, = 0, 



(E2) 



where 



/3(x) = V • €, 



(E3) 



the radial component of £ is £ r e r , and the functions gi(r) = Y\pjp 2 and g%{r) = Y\pA r jp are fixed functions 
determined by the background stellar model and by the adiabatic index Ti of the perturbations. The last two terms 
involving the function vanish for zero-buoyancy stars. 

The argument of Lockitch and Friedman |57|, simplified for Newtonian stars, is as follows. The perturbed pressure 
p = p + Sp and perturbed density p = p + Sp obey the equation of hydrostatic balance and Poisson's equation, to 
linear order in the perturbation. Hence, the pair of functions {p, p) describe a static self-gravitating perfect fluid 
configuration, if one assumes that to each solution of the linearized fluid equations there corresponds a solution of 
the exact equations. However, it is known that all static, self-gravitating perfect fluid configurations are spherically 
symmetrical Hence, all perturbations (Sp,5p) must correspond to displacements to nearby static spherical config- 
urations. Solutions of this type include (i) perturbations taking the star to a nearby static, spherically symmetric 
configuration with a different total mass M. If one assumes that the change SM to the total mass of the star vanishes, 
such solutions are disallowed, (ii) Solutions with I = 1 of the form £ = const corresponding to the displacement of 
the center of mass of the star. If one demands that the perturbation leaves fixed the center of mass of the star, these 



26 Note that the spherical-symmetry theorem does not require the barotropic assumption VpxVp = 0. It requires only that 
the perturbed pressure p and perturbed density p satisfy the equation of hydrostatic balance and Poisson's equation. 
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solutions are disallowed, (iii) Solutions with I = where the background star is marginally stable to radial collapse. 
Henceforth we will assume the background star is such that there are no I = zero-frequency solutions of this type. 
If one makes sufficient assumptions to outlaw the cases (i), (ii) and (iii) above, it follows that all solutions to Eq. 



(E2) must satisfy 



Sp = Sp = 0. 

We now consider the cases of zero-buoyancy and nonzero-buoyancy stars. 



(E4) 



1. Zero-buoyancy stars 

We showed above that all solutions in Ho have Sp = 0. However, for A = 0, Eq. ( fE2| ) is always satisfied when 
Sp = 0, since the last two terms vanish. Hence, the space Ho consists precisely of the perturbations with 



Sp = -v • (pO = o. 



(E5) 



Note that this characterization does not hold for rotating stars, where there are zero frequency modes with Sp =/= 
(see Appendix ^) . 



We next derive a property of Ho that is used in Sec. [II above. Let Po to be the operator that projects orthogonally 
onto the subspace Ho- We now show that 



P • | = if and only if V x £ = 0. 



(E6) 



Let T lJ (x) be any antisymmetric tensor that vanishes on the boundary of the star. Then an integration by parts 
shows that 



/ 



d 3 xVuaT^ = (^,v) 



(E7) 



where the vector v is given by 



P 



(E8) 



Now if Pp ■ £ = 0, then £ is orthogonal to all vectors v in Ho- But vectors of the form flEq) are automatically in Ho, 



by Eq. (E5). Hence both sides of Eq. (E7) vanish. Since this is true for all T lJ (x) that vanish on the boundary, and 
since Vj£j is continuous, it follows that V x £ = 0. It can be checked that the argument also carries through in the 
other direction. 



2. Non-zero buoyancy stars 



Turn now to stars with buoyancy forces. For these stars, the Lockitch- Friedman argument shows that the space of 
zero frequency modes is prec isely the set of axial vectors (except possibly in the / = sector). 

We use the expansion (ful) of La grangian perturbation on a basis of vector spherical harmonics. Here the first two 
terms are the polar part of £ and the last term is the axial part. It is easy to see that th e ax ial part gives a vanishing 
contribution to both Sp and = V • £, and hence all axial vectors lie in Ho, from Eq. (E2). Thus the C; m terms do 
not contribute at all to Eq. (E2) and so we can drop them and focus on the A\ m and Bi m terms. The vanishing of 



the density perturbation Sp — —V • (p£) from Eq. (E4) can be written as 



-(pA lm y + 2A lm /r - 1(1 + l)B lm /r 2 = 0. 



(E9) 



Also, the relation between the Eulerian density and pressure perturbations for A ^ is, from Eq. (2.4), 



8p=ElP(5p + pH-A) 
P 



(E10) 



Substituting from Eq . (E4) now shows that the radial component of £ vanishes. Hence Ai m = for all l,m, and it 
follows from Eq. (E9) that Bi m = for all l,m with I > too. Since we are restricting attention to stars for which 
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there are no I — zero- frequency modes (see above) , it follows that the space Tio of zero frequency modes is therefore 
precisely the set of axial vectors. 

Lastly, let Pp b e the orthogonal projection operator that projects onto Ho- In this case it follows from the definition 
of axial [Eq. (|5.l| ) with Ai rn = Bi rn = 0] that 

P • £ = if and only if r • (V X £) = 0. (Ell) 



APPENDIX F: NON-EXISTENCE OF MODES WHOSE FREQUENCIES SCALE AS THE SQUARE 

ROOT OF THE ROTATIONAL FREQUENCY 



In this appendix, we show that there are no modes whose frequencies and mode functions can be expressed as 
power series in the square root \/f2 of the star's angular velocity and whose frequency in the non-rotating limit is 
zero. As discussed in Sec. [II A above, one might suspect the existence of such modes from the general perturbation 
theory analysis of [M and the fact that non-rotating stars have Jordan chains of length one. The result is valid for 
both zero-buoya ncy and nonzero-buoyancy stars. The non-existence of such modes is related to the vanishing of the 
matrix element (D7) discussed in Appendix |d| above. 

We assume a one parameter family of modes £(fi),u/(n) with expansions of the form 



and 



1(0) = £(°) + VfU (1/2) + + o(n 



3/2\ 



(Fl) 



(F2) 



If we substitute the expansions (|Fl|) -(|F2|) and also the expansions ( 3.2 )— ( 3.3) for the operators B and C into the 
quadratic eigenvalue equation ( |2.18| ) we get at orders O(O ), O(VU), and O(fi) the equations C' ' • £ (0) = 0, 
C(°) • = 0, and 



C (0) - W ( 1 /2)2|(0) =Q 



(F3) 



We now multiply Eq. ( |F3| ) by the projection operator Pp into the space Ho of zero-frequency modes, and use the fact 
that P • £ (0) = £ (0) and P • C^ ) = 0. This gives 



w 



(1/2) 2 AO) 



0, 



(F4) 



which shows that u/ 1 / 2 ) = 0. 



APPENDIX G: DERIVATION OF SELECTION RULES 



In this appendix we derive the two axial selection rules (5.14) and (5.15) and the selection rule (5.16) for an 
incompressible star. Throughout this appendix we will drop all terms that are not of zeroth order in £1. Therefore we 
will use the pressure, density and potential profiles of a spherical background star. 



1. Axial selection rules 



We assume that the modes £ B and £ c are axial, and therefore 

V • $ B = V • $ c = 5 PB = Sp c = & = & = 0- (Gl) 

In addition we assume that the Eulerian density perturbation 5pA of mode $, A vanishes, and so from Eq. ( |2.3|) we 
have 



(G2) 
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Inserting these results into the expression (4.2C) for the coupling coefficient we obtain 



2kabc = Jd 3 xp [(Ti - 1)(V ■ £ a )Zbc + Xabc + Xacb] - J d s x P C A ? B i 



1 c3 ck c^k 



Now the background potential is a function of r only, so we have 



1 d (1 



>;jki 



dr V r dr 



(xiSjk + x k 5 lJ + xjSik) + 



1 d 

r dr 



1 d f 1 dcj) 
r dr \r dr 



XiXjXfc. 



(G3) 



(G4) 



Substituting Eq. (G4) into Eq. (G3) and using Eq. (Gl) gives 

2kabc = J d 3 x p (xabc + Xacb) + J d 3 x p(T x - 1)(V • £ A )E 



BC 



d i p- 



d (\ 



dr V r dr 



(G5) 



We now evaluate the first term in Eq. (G5) using the definition (4.25) of xabc- This gives 
Jd 3 xp X ABC = Jd 3 x i>Cu XU : iXr : , J d 3 xp{C A e B ;k)dc ^ 



(G6) 



where the second form of the integral was obtained using Cs -fe? = ^B-jk ~ 0- Upon integrating by partsQjwe find 



J d 3 X PXABC jd 3 .,- fa&tfg ;i + P&& :ki 



k \ 
C;ij> 



We can rewrite this using the relationship that follows from p = p(r) 



the result is 



d x pxabc = d x 



PtikCA^B^C ;i ~ pCa£b ;k^C ;i 



(G7) 



(G8) 



d^a; 



(G9) 



wher e w e have again integrated by parts. We can evaluate the term involving p-jki using an equation analogous to 
Eq. (G4), which gives 



d 3 x pxabc =~ d 3 x 



dr \r dr 

Next, we symmetrize over the indices B and C. This gives 



X (iX ) & (O ' €b) +P-,Ma-AUc +P-MaZ B Jc +P&&;k&fi 



(G10) 



J d 3 x p (xabc + Xacb) = -2 J d 3 x 



d fldp\ 
dr \r dr J 



Ca^a^b)+P-,Ma-AUc 



d 3 x 



P^B 4C + P&& ;dc;ij + P'M^C + P&& ;£b 



//•' 



(Gil) 



27 Throughout this appendix we discard surface terms that arise from integrations by parts which involve the radial component 
£^4 of the mode function £ A evaluated on the star's surface, because the quantity 



Ap 
P 



Tip dp 
cr dr 



Ti 



4>'& 



must vanish on the surface fkt 5| 
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where on the right hand side we have relabeled some indices. We can rewrite the second line of Eq. ( |G11 ) as 



k 



d 3 X 



P-JkiCv^Uc + P:jkC A AUc + /'• jXc,, + pCa An jAr- 



CO c k 



(G12) 



where we have agai n inte grated by parts. Rewriting the terms invol ving p-jki and p-jk the same way as before, and 
using the definition ( 4.24 ) of ~E.bc, we can write the expression ( Gil ) as 



d 3 x p(xabc + . 







_dr 


\r dr J 



,c 



+ J d 3 x [ P 4 A z BC +pe A .^Bc]- 

Now since the Eulerian density perturbation — (p£ A );i vanishes, we have 



dp \ amp 



-p(T-l)&, 4 , 



(G13) 



(G14) 



where T is the adiabatic index of the background spherical star. Combining Eqs. (G5), ( G12| ), and (G13) now yields 
2^ abc = J d 3 x (€ B -€ C 

+ J d 3 xp(T 1 - T) (V • £ A ) E BC 



(G15) 



Using Eq. (G2) and the eq uatio n of hydrostatic equilibrium in the background star one can show that the first line of 
the right hand side of Eq. ( p!5| ) vanishes. Therefore we find 



kabc = \J d 3 xp{T 1 - T) (V ■ £ A ) Zbc- 



(G16) 



We can derive two different selection r ules f rom Eq. ( G16| ). First, if £ A is axial, then V ■ £ A = and therefore 
kabc = 0; this proves the selection rule (5.14) for three axial modes. Second, if the perturbations obey the same 
equation of state as the background star, then the adiabatic index Ti of the perturbations coincides with the adiabatic 
index T of the background star, and k A bc again vanishes. This proves the second axial selection rule ( p. 15 ). 



2. Incompressible selection rule 

We can co nside r an i ncom pressible star to be the limit T\ — > oo of a zero-buoyancy polytrope with T = Tx. The 
expressions (4.20) and (4.27) for the coupling coefficient kabc are valid for incompressible stars, since they are valid 
for finite Ti and they have a finite limit as Ti ^ oo. Note however that one must be careful to include all the 
distributional contributions in the integrands at the stellar surface that arise in the Ti — > oo limit. 

Consider now a linear mode £(x) of the large-I^i star, whose ri — * oo limit is a mode of the incompressible star. 
Clearly the Lagrangian perturbation of the density must scale as 



which implies from Eq. (114) that 



Ap l_ 

p 00 IV 



V-Coc— . 
1 1 



It follows that to zeroth order in \/T\, the Eulerian density perturbation is, from Eq. ( [2.3] ) 

*P = -V-(p£) = -€-Vp-p(V-e) 

= Cf>S(r - R) 



(G17) 



(G18) 



(G19) 
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where p = 3M/47ri? 3 and M and R are the stellar mass and radius. 

Consider now three modes £ A , £ B , £ c that satisfy the hypothesis of the selection rule ( 5.16|) . That is, th ey all have 
vanishing Eulerian density and pressure perturbations to zeroth order in 1/Ti. It follows from Eq. (G19) that they 
all have vanishing radial components at the stellar surface, 

C(r = R)=0. (G20) 



Consider now the expression (4.27) for the coupling coefficient. In this expression, the second, third and fifth terms 
vanish since V • £ = for each mode, and the fourth and sixth terms vanish since each mode has vanishing Sp, S<j) 
and Sp. The remaining terms yield the following expression for the coupling coefficient: 

KABC = -~ / d 3 x£ A £ 3 B & \P;ijk + P<j>;ijk] • (G21) 



The integrand in Eq. (G21) vanishes in the interior of the star, since p and </> are quadratic functions of position for 
a constant density star. However, there are distributional contributions to the integral at the stellar surface, since p-i 
and (j>.ij are discontinuous across the surface. To evaluate these contributions, substitute in the equation of hydrostatic 
equilibrium in the background star, + pcj)-i — 0, which yields 

kabc \J d?xC (A e B ^C) PPiritik + <t>;iP;jk] ■ (G22) 

The first term vanishes since ip.ji is finite on the stellar surface, and for any of the modes we have 

e*Pifc = ~H r (r =R)=0 (G23) 
from Eq. ( G2C| ) . Integrating the second term by parts gives 



KABC = / 



P-.k- (G24) 

i 



When the derivative of the square bracket is expanded out, each term will contain either a factor Cp-.i or a factor 
both of which vanish on the stellar surface. Therefore we obtain 

kabc = (G25) 

to zeroth order in 1/ri and in fl. Note that the cancellation of the distributional components in the expression ( |G2l| ) 
would not be obtained if one uses the Cowling approximation of neglecting the gravitational terms in the coupling 
coefficient, as in Ref. p8|. 



APPENDIX H: CALCULATING COVARIANT DERIVATIVES USING THE SPIN WEIGHTED 
SPHERICAL HARMONICS AND THE 3 OPERATOR 



We define the "spin- weight" s of any tensor as in |g7j (also see |36[ for a very readable refe rence) , in terms of how 
it transforms under rotations around the radial direction. For example, the vector m of Sec. VI B transforms under 
these types of rotations as m — > m', where 



m' = e^m 



where ip is the angle of rotation. We say that this vector has spin-weight s 
have spin-weight s if it transforms as $ — > *&' where 



(HI) 

1. In general, a quantity $ is said to 



(H2) 



In particular, the basis vectors l,m and m defined in Sec. [VI B| have spin- weights 0, 1 and —1 respectively. Thus we 
can decompose any given vector function, B(r,0, if), into three components with definite spin-weight by dotting it 
with one of the three basis vectors (1, m, m). Each scalar component can then be expanded in the complete set of 
basis functions s Yi m , e.g. 
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B-l = CoM,^)=X>oM o^a 

A 

B- m = ( 1 (r,e,^) = Y / bi(r) iY A 

A 

A 

For notational brevity we note that m and m, and consequently £1 and are complex conjugates of each other. 
Thus, in subsequent calculations we will write +c.c. when convenient.^] The connection coefficients on this basis are 

pi _ _ pm _ 

mm Im 

r m ^rm ]_£Ot0 

mm mm \/2 ^ ' V / 

where all the other components are zero except for complex conjugates of the above, e.g. 

/rim \* pm _ co ^ ^ fHSi 

\ mm) mm \/2 ^* ' V / 

It follows that the components of the covariant derivative of B are 



B> 





= B l , 




B -m 


= B m m + 


1 , cot 6 
-B l + —=- 

r s/2r 


a ;m 


r>m 

— B ,m ~ 






= B m , 




B ;m 


= B jn ~ 


-B"\ 

r 



(H6) 

where all others are zero except for complex conjugates of the above. Now consider the 3 ("edth") and 3 operators 
of Ref. p7|], defined by their action on a function, £ s of spin- weight s, 



Q( s = -(sine»)V2(m- L)(sin6l)~ s C s = - (d e +icsc9d v - scot0)Cs 

ECs = (sm60~ s V^(m-L)(sin<9) s G = (-d e + % esc 0d v - s cot 0) (H7) 
Using these definitions in Eq. flH6| ) we get 



1 1 

\J2r r 



and complex conjugates of the above. The 3 (3) operator, when acting on the spin-weighted spherical harmonics, 
raises (lowers) the spin-weight by one, 



,Y lm = y/(l-s)(l + S + l) s+1 Yi m 

sYlm = -y/{l + 8)(l-8 + l).-xYl m . (H8) 



28 Note that in Section VI B we choose to implement this by simply adding terms with s replaced by — s. 
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Thus when calculating a covariant derivative such as B m m we can write 



(H9) 



and similarly for the other covar iant derivatives. This allows us to define the functionals as in Eq. (3.15). For 
decomposing vectors of the form (3.1), note that 



m-V Y A 
m-V Y A 



A 

x/2> 1 
A 

A 



1a 
iYa 



m • L Ya = -p \Yk 

m L Ya = --j= -\Y A . 



(H10) 



Using these relations together with Eqs. (5.6) and (5.7) we find the functions /q , and /_ for z— parity odd modes 



/TO 



/TO 



/TO 



Jj+m+l,l 



Vj+ m +i + m + l){j + m + 2) 

Uj+ m + m)(j + m + 1) 
r 



Oj+m+1,2 



V j+m+ i + m + l)(j + m + 2) 
r v 7 ^ 



Uj+m + + TO + 1) 

r V2 



and for z— parity even modes 

/TO 
f+ l (r) 

/TO = 



W j+r 



Vj+ m a/(j +m)(j + m+ 1) 
r y/2 



5 iH 



+ m + l)(j + m + 2) 

i ;= Oj+m+l,J 



72 

U/+m \Z(j + + m+ 1) 
r V2 



-\-mJ 



Uj+m+i ^/(j + m+l)(j + m + 2) 
r V2 



(Hll) 



(H12) 



APPENDIX I: DIRECT DERIVATION OF SECOND ORDER EQUATIONS OF MOTION 

There are at least three ways to derive the equation of motion for the second order Lagrangian perturbations: (1) 
Expand the variational principle to third order in the Lagrangian displacement £, and then vary with respect to £; (2) 
Vary with respect to £ to get the equation of motion, then expand to second order; (3) Directly perturb the equation 
of motion to second order in £: 



(II) 
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The second method is carried out in Appendix Here we carry out the third method. 
Recall that the Eulerian and Lagrangian perturbations of a quantity Q are defined as 

5Q = Q(x J *)-Q (x,t), (12) 
AQ = Q(x + £(x,i),t)-Q (x,t), (13) 

where Qo is the value of Q in the unperturbed solution. [In this appendix we shall sometimes use po, pa and (/>o to 
denote the background density pressure and potential; these quantities were denoted simply p, p and <f> in the body 
of the paper.] Suppressing the time dependence of the quantities, we can write to second order 

AQ = Q(x + 0-Qo(x) 

= Q(x + - Q (x + + Qo(x + £) - Qo(x) 
= <5Q(x + + Qo(x + 0-Qo(x) 

- <5Q + fV^Q + fV.Q + ^e^ViVjQ. (14) 

In the last line we have dropped the subscript on the last two terms. Note that SQ in the second term on the 
right-hand side need only be evaluated to first order. Substituting the first order relation 

SQ = AQ - CViQ (15) 

in this term gives the alternative relation 



AQ = SQ + rv 4 AQ + C V»Q - CV,£ J V,Q - V l V ] Q. (16) 



By definition of the Lagrangian displacement vector the unperturbed and perturbed fluid trajectories are related 

by 

4 = ^*), a* =<:<(*) +£V(t),t). (17) 

Thus 
Hence, 

A„. = f , (19) 

an equation that is exact to all orders in £ (see also Ref. (42|). 

An expression for Ap follows from conservation of mass. The variation in the total mass due to the perturbation is 



= / p(-x)d 3 x- / p (x)d 3 a; 

JV+AV Jv 

= f [p(x + £)J-po(x)]d 3 x. (110) 
Jv 



Note that we have defined the Lagrangian perturbation of the vector v % by subtracting the components at two locations 
separated by f . This procedure is valid for Cartesian vector components in flat spacetime. We can derive formulas valid in 
arbitrary coordinate systems by working first in Cartesian coordinates, and then rewriting all partial derivatives as covariant 
derivatives at the end. 

Friedman and Schutz have suggested defining the Lagrangian perturbation by mapping the perturbed vector back to the 
original location so that vectors are subtracted only at the same point (Lie derivative) . This gives different intermediate values 
for the Lagrangian perturbations of various quantities than the procedure followed here. However, when the final equations of 
motion are written in terms of £, the results are the same since one is using the same definition of £. While the Lie derivative 
is essential in relativistic applications, the present approach is simpler here. 
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Here we have made the transformation x' = x — £(x) in the first integral to make the volumes the same, and then 
dropped the primes. The Jacobian of the transformation is 



J = ^xT = det ^ + L ^ ) ' 



where 



3 ~ dxi ■ 



Setting 5M = and substituting p(x + £) = pa(x) + Ap in Eq. (110), we get 

Ap I- J 



P 



J ' 



exact to all orders in £. Keeping terms through second order from Eq. (J19), we find 

— = -Vif + kvif ) 2 + V^ fc V fc f ]. 
P 2 



(111) 



(112) 



(113) 



(114) 



Consider an equation of state p = p{p, p) where p denotes quantities such as entropy or composition. For adiabatic 
perturbations in which the composition does not have time to change Ap = 0, so to second order 



op 



1 d 2 p 
iu 2 °P 



In terms of the adiabatic index governing the perturbations 

dhip 



d\np 



we get 



Ap Ap 1 



Ap 
P 



(115) 



(116) 



(117) 



Here Ap need be evaluated only to first order in the last term. 

To derive the equation of motion for the perturbations, we need to perturb various derivatives of quantities. Consider 
for example dQ/dt. By Eq. 



dQ s dQ 



: dQ 



dQ 



dQ 



at at * dt * ot ' 2~ 

The Eulerian perturbation operator commutes with partial differentiation: 

* at = dt 6 ' 



dt 



so the first term on the right-hand side of Eq. (118) is 



d 



d 



dt ^ dt 
where we have used Eq. (M) again. Thus 



— SQ = — [ AQ - CVt&Q - fViQ - ViV^Q ) , 



ae 
at" 



(ViJQ + ViQ + ^'ViVj-Q) 



(118) 



(119) 



(120) 



(121) 



Here we have used the symmetry of VjVj to rewrite the last term without the factor of 1/2. An alternative form of 
Eq. (121) follows from substituting Eq. (M): 
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dO d dP 
A-^ = —AQ ~ JL (V,AQ + V t Q ~ V^VjQ) 



Of 



(122) 



A similar equation holds if d/dt is replaced by d/dx 1 , and then we can let d/dx l — » to get a tensor equation valid 
in any coordinate system. 

For the perturbation of the acceleration, we have 



A^ = A 
dt 



d_ 

dt 



u k V k ) u t 



A— u r + u k AV k Ui + Au k V k Ui + Au k AV k Ui 
at 



(123) 



Note that our formulas hold in non-Cartesian coordinates whether we perturb the covariant components Ui or the 
contravariant compo nen ts u\ for the reason explained in footnote ^9|. Commute the operators Ad/dt and AV k in 
Eq. (123) using Eq. (121) and the corresponding equation for Vfc. Many terms cancel, and we are left with 



, duj d 2 P 
A — - = — - 

dt dt 2 

For a time-independent unperturbed star du % /dt = 0, and so 



dui 
^~dt 



d_ 

di 



d_ 
m 



u k v k e 



d 2 f dP 



Using Eq. ( |l4"| ) and commuting 5 and Vj, we get 

AV.,0 = V 4 <5</> + ^VjVi5(j> + £%-V< 

The Eulerian perturbation is found from Poisson's equation: 

V 2 = AnGp =>■ V 2 

which implies 

Spf 



= 4nGSp, 



-G 



d 3 x'. 



We get an expression for Sp by substituting Eq. (114) in Eq. (|^): 



After some algebra, we can simplify this to 



(124) 



(125) 



(126) 



(127) 



(128) 



(129) 



^=-V l ( P e) + -V 4 V J (pfe)- 



(130) 



We have 



30 Note that the second term in Eq. (130) typically includes a 8- function at the stellar surface coming from the sec ond derivative 
of the density. This 8- function is physical and should be included; it corrects for the fact that the expression (130) vanishes 
outside the surface of the unperturbed star whereas in reality Sp can be nonvanishing there. 
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1. 



1 



A - Vip = A - Vip + -AV lP + A - AV ( p 



P 

-\a p +\{A p ) 2 
p 2 - p A 



ViP + -AV lP - X^ApAWip. 
p p z 



(131) 



Substitute Eq. (122) with V^p replacing dQ/dt: 



pA ( - ViP 



V# + V. ( Ap 



V^V.-Ap V,,rT,/' • V,,,: 
__^( Vi Ap- Vi^'Vj-p) 



(132) 



Recall the definitions of the tensors J 'j and 3 J 'j and their traces from Eq. (fh9). Their divergences satisfy the identity 



vM-v^e^-Vits-e). 



(133) 



Substitute Eqs. (|I14j) and QI17D in Eq. (|I32| ) to get 

/'A (Vp) - Vjf - Vif'VjP - V,;/,l ,V,e':' 



+-(e-s)ViP 



i 



i 



pri-(e + H) + -pG rf-r x +p 



dp 



+(s fc , - e fc i)v fe p + v^v^prxVfe^) - v,'T,i,»iT t .f) 

It will be convenient to rewrite the second order terms on the right-hand side of this equation as a divergence. (The 
first order terms can also be written as a divergence.) The last two terms are 

Vi^Vj-CpiWfcf*) - WVi(priV fc £ fc ) = v^vrfViV*^) 

-v^v^piw*^) 
= Vjipr^i) - Vi(pr x e). (i34) 



Also, by Eq. (133) we have 



Thus Eq. ([1341) becomes 



(s fc , - e\)v kP = v fc [p(s fc i - -pv fc (s*< - e fe 
= v fc [p(s fc J -e fc 4 )]-ipv l (s-e). 



paQv.pI : V ( i'V,p- Vit-'V,-p- VMyTjV.i'i 



(135) 



+v 



^p^-^i + ^ip^-lf+p 



dp 



e 



+pH^ + ^<p' iP (ri - i)h}. 



(136) 



Multiplying Eq. ( |ll| ) by p, and using eqs. (124), (|36j), and (|3C|), we obtain the equation of motion for the pertur- 
bation to second order: 

= -Vi^'ViP + VifVj-p + V,i/-I ,V ; i-'i - pV^^V - pf'VjVi^ 



-VijpfT! - 1)0^ + i^ipjCT! - l) 2 + P ^ 

+ps j , + ~*\p(ri - i)s} - pv 



e 



(137) 



50 



Here we have split the quantity 5(f> defined by eqs. (128) and (130) into first and second order pieces 8W<j> and 5^(j). 
Accordingly the first line on the right-hand side of Eq. (137) gives the first order terms in the equation of motion, 
while the remaining lines give the second order terms. 

Finally, we note that the above derivation was carried out in the incrtial frame. To obtain the equation of motion 
in the rotating frame, we replace the inertial-frame Euler equation ( |Tl| ) at the start of the calculation by the rotating- 
frame Euler equation (2.2). This involves adding to the left hand side the term 



and replacing the background po tent ial <f> by 
Taking the variation of the term (138) using Eq. 



2S1 x u, (138) 

^rotj where the centrifugal potential (f> m t is given by Eq. (pT0| ). 
yields 



2n x 



dt' 



(139) 



which should be added to the left hand side of Eq. ( 137 ) , and substituting 0^0 + (j) Iot on the right hand side 
results in the additional term — p£ J Vj Vi</> ro t on the right hand side. In the application of this paper, the background 
velocity in the rotating f ram e vanishes, and so we can replace d£/dt with d£/dt and d 2 £/dt 2 with d 2 £ /dt 2 . With 
t hes e mod ifications, Eq. (137) reduces to the result quoted in the body of the paper, given by Eqs. (|1|), (|2.9| ), (|f.l|), 
(O), Q and (O). 



APPENDIX J: VARIATIONAL DERIVATION OF SECOND ORDER EQUATIONS OF MOTION 



Derivations of the hydrodynamical equations directly from a variational principle have a long history |67| |70| , and 
have been applied to perturbation theory for nonrotating stars (e.g. Ref. @). Hydrodynamical equations can be 
derived from a Lagrangian viewpoint by extremizing the action 

S = J dtC{t)= J dt[T{t)-V{t)], (Jl) 

where the Lagrangian is the difference between kinetic and potential energies for the fluid. Included in the potential 
energy are both the gravitational potential (both the self potential and external potentials) and the thermal energy. 
The specific variational principle we use here is the following 69 1 . Label fluid elements by their positions yo at some 
initial time t — to, and let the position of fluid element yo at time t be y(yo,i)- Let po(yo) be the mass density at 
time to, then the element of mass in the fluid is dM = pod 3 yo — pd 3 y where p = pa/J and J is the Jacobian 



Then the kinetic energy is 



' 'II dM (^) M-tl 

yo 



2 / v at 



and the potential is 



r - I dM£(p) - ° / *» J ^(vMM [ MM ■ (J4) 



Here p is understood to be the functional of y(yo, t) given by p — po/ J, and £ is the internal energy per unit mass of 
the fluid, related to the equation of state by dS/dp — p(p)/p 2 . With these definitions, varying the action (|jl| ) with 
respect to y(yo,t) gives the Euler equation p9| . 

We can generalize the above variational principle slightly by allowing the energy per unit mass to depend on 
additional hydrodynamic variables such as specific entropy or composition. We denote the set of such variables by 
the vector p, so that 

£ = £{ P ,p). (J5) 

In this case, the action principle ( |jl| ) does not determine the evolution of the variables p. However, if we restrict 
attention to situations in which the variables pi are determined by the function y(yo,£), then the action jjl| ) does 
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determine the motion of the fluid. An example of such a situation is where each fluid element conserves its values of 
ft; see below. 

We now apply this principle to perturbations of a uniformly rotating star. We denote the background quantities by 
Po, pa and <j)Q] these were denote by p, p and (f> in the body of the paper. We can reformulate the variational principle 
( |jl| ) by transforming to the co-rotating frame. Define Vp = U(t) • xo and y = U(i) • x, where U(t) is the appropriate 
orthogonal rotation matrix. The form of the piece ([jj) of the action does not change under this transformation, so 
we merely have to transform the kinetic energy (J3). In the rotating frame, the position at time t of a mass element 
with comoving coordinate xq is 



x = x(x , t) = x + £(x , t). 



(J6) 



[Note that the xo of this appendix was denoted x in the body of the paper.] The position as seen in the inertial frame 
in which the star is at rest is 



y% = Uij(t)[x 0j + £j(x ,i)]. 
The inertial-frame fluid velocity is found by differentiating this equation with respect to time: 

••■-^[xoj+^,t)] + U v (t)^ M 



dt 

Using the identity dUij/dt = ei P q£l p U q j(t), we find 

V- 1 {t)-u = nx [x + £(x ,t)] 
where we have used Q m = U pm (t)Q p . This yields 



at 



^(XQ,t) 
dt • 



|u| 2 = |Ox [x +£(x ,i)]| 2 -2[x + £(x ,<)] • Ox 



ag(xo,t) 

dt 







)+ 


dt 



(37) 



(38) 



(39) 



(J10) 



To find the equations of motion for the perturbations of the star, we can dispense with all terms in the Lagrangian 
that are either zeroth or first order in £ and its time derivatives^ We will denote with subscripts 2 quantities from 
which the zeroth and first order pieces in £ have been subtracted. Thus, the relevant piece of the kinetic energy is 



To 



dM 



dt 



2^ 2 |£l 



d£(xp,f) 
dt 



We also note that we can define the canonical momentum 

dT 2 d£ 



d(d£/dt) dt 



where 7~2 is the integrand in Eq. (Jll), and that 



A Hamiltonian for the system is thus given by 



H 9 



dM 



■K-(flx£) 



\n- n x £1 



|7r-fJx£| 2 \Slx£\ 



\nx£\ 



v 2 (t), 



(Jll) 



(J12) 



(J13) 



(J14) 



where V% is the potential (J4) with zeroth order and first order terms dropped, which only depends on £, Vo£ and 
the background star, but not on 7r. Hamilton's equations then yield the equations of motion 



31 Note that there is a term proportional to j in |u| 2 , but varying the action with respect to this term does not contribute to 
the equations of motion. This is because the result of varying this term is a total time derivative. 
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£ = 7T - ft X £ 

f ~ 

7T = It X 7T, 



(J15) 



where the accelerations f is determined from the variation of the potential V 2 - In the linear approximation, these 
results are equivalent to the first six equations in Appendix A, with no external force. We can combine Eqs. (J15) 
into the single equation 



<9 2 £ d£ f 

« +2 nx-| = --nx nxfl. 

at 1 at po 



(J16) 



It remains to compute the variations of the potential, V2, which contains two parts, internal energy and gravitational 
energy. The internal energy is 



U= / dM£{v,p) 



(J17) 



where v is the volume per unit mass. (We express £ here as a function of v instead of p = As explained above, 

p is a set of variables such as the specific entropy and composition of the fluid. The volume per unit mass in the fluid 
is related to its value in the unperturbed state, i>o, by 



J(x,x )z; , 



(J18) 



where J(x,Xo) = det(l + L), with = <5jj + d^/dxoj [cf. Eqs. (Ill) and (112)]. The Jacobian J is given by the 
exact formula 



J(x, x„) = 1 + V • € + i (0 - 5) + i(V ■ £) (0 - 35) + 



(J19) 



where 0, 5 and \ & r e the traces of the tensors defined in Eq. (4.6). This formula is essentially just the characteristic 
equation for the 3x3 matrix <f' ._■ . 

To go further, we need to specify how p changes in the presence of fluid displacements. We assume that for each 
fluid element p remains frozen at its background value po [3 This means that we focus on perturbations that preserve 
for each fluid element the composition, specific entropy, etc. of the unperturbed star. We are essentially assuming 
that the timescales for the fluid to relax to local composition are relatively long] 33 ] 

Although it is easy to subtract the contributions to U that are zeroth and first order in £, for computing its variation 
it is easier to work with U rather than U2, and subtract off the displacement field through J(x, Xq) only. The result is 



5U = - J dM v p[vo^(x,x ),/io]6J(x,xo) 

d 3 x p[u J(x,x ),^o]<W(x,xo), (J20) 

wherep(v) = —d£(v, p)/dv is the pressure. To simplify further, we substitute <5J(x, xo) = J(x, XQ)(l+L)~ j i dS^i/dxoj, 
and then integrate by parts, and impose the surface boundary condition p = 0; after subtracting the pressure gradient 
force in the background star, the result is 



6U 2 = / d 3 x 6k 



d 



dx . 



p[v J(x,x ),/loR(x,xo)(l + L)^ - %p(w ,^ ) 



(J21) 



Note the exact relation 



J(l + L)£ = [1 + V • € + i(0 - 5)]^ - + Eji - Oji. 



(J22) 



32 More generally, we could assume that p changes only in response to gradients in £, and, since the various components of p 
are scalars, presumably would depend only on scalar quantities that can be constructed from Vo£. 
33 To do better, we could introduce equations of heat transport, composition and so on. 
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(For a derivation, see e.g. Ref. jn].) Expanding the quantity in square brackets in Eq. (J21) to second order, we find 



p[uoJ(x,x ),//o]J(x,x )(l + L) Si - 5 i3 p(v , jl ) =Po(^o,Mo)i ^(1 - Ti)Vo • £ 



dx 0i 



2 



ehr, -Mr, - lie 



(J23) 



which gives the same result, to second order, as Eq. (137). Here, we have defined a generalized version of the adiabatic 
index of the perturbation, 



Ti = - 



d In po (vo, po) 



d\np (p ,Po) 



(J24) 



OIiivq J ,-r V 91np „ 

which also applies in a zero entropy star, for perturbations whose composition, and possibly other characteristics, 
remain fixed at their unperturbed values. It is straightforward to extend the expansion to any desired order, and also 
to expand the internal energy per unit mass, £(vqJ, /2o) , to any desired order. 
Finally, we need to consider the gravitational potential energy, 

dMdM ' (J25) 



2] |x -x' +£(x ,t)-£(x(„f)| 



Varying with respect to the displacement field gives 

dV G = 17 / dMdM d£(x ,i) - 0£(x ,i) — 

27 l( x o -x + £(x ,<) -|(x ,i)| 3 

by switching primed and unprimed variables, and subtracting off the zeroth order gravitational force, we find 

x -Xo + £(x ,£) -£(xo,t) x -x 



5V G2 = G J dM6£(x ,t) ■ J AM' 



(x -x^+^(x ,t)-€(x^i) 



l x o- x ol 3 



(J26) 



(J27) 



Gathering results, we see that the acceleration in Eq. flJlq ), computed to any desired order in the amplitude of the 
Lagrangian perturbation field, £(xq, t), is f/po = fp/po + fa/po, where 



f j = ^_ 

P dx 



Po 



:j l 
G I AM 



p[uoJ(x,x ),^o]^(x,xo)(l + L)jl - <%p(vo>#o) 
xo-Xo+C(xo,t)-^(xQ,t) x -x 



|(x -x + £(x ,0-£(x ,i)| 3 |x -x 



' 13 



(J28) 
(J29) 



Equations (IJ23J) and (J28) allow us to calculate fp to second order in the perturbation; the result is Eq. (k 



Simi larly , let us expand = 



(i) 



F (2) 



Eq. (J2E). To first order, the result of the calculation is that 



the contributions from various orders can be found by expanding 



-Po£ 7 (xo,O v o,iV Oj 0o(xo) - p o V o ,i(5 (1) 0(xo,i), 



where 



5 (1 V(x ,t) =G / d 3 x 



, V > (x )^(x ,i)] 



x 



The second order contribution to the gravitational force is (compare Eq. ||I37[ ) 



where 



,(2) 



5^<t> : 



poijiu 



G 
2~ 



/ d 3 x 'p (x )^.(x , 0) t)^(x^,t)V OJ V Oifc |x -x 



or 1 - 



(J30) 
(J31) 

(J32) 
(J33) 



It is possible to integrate by parts to change the form of 5^4> to show it coincides with the expression (4.4) 



An integration by parts is required to get this form of the gravitational force. The outer boundary has been taken to be 
outside the star, so if the mass density goes to zero discontinuously there, a surface term must be included. 
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APPENDIX K: ENERGY AND ANGULAR MOMENTUM OF STELLAR PERTURBATIONS 



In Sec. K 1 of this appendix we give general formulae for the physical energy E p h ys and physical (z-componcnt 
of) angular momentum J p h ys of a perturbation, and also the rotating-frame canonical energy -E C an,rot ^fj- General 
expressions for these quantities have been derived by Friedman and Schutz pl| , in a context more general than that 
considered here (uniform rotation of the background star). However, the results of Friedman and Schutz are expressed 
in terms of an inertial-frame Lagrangian displacement £ in . For a uniformly rotating background star, it is more natural 
and simp ler to use instead the rotating-frame Lagrangian displacement £ as we do in this paper. The expressions in 



Sec. K 1 below can be obtained by specializing and translating the results of Friedman and Schutz, or, as we outline 



below, can also be easily obta ined directly. 



We also derive, in Sec. |K 4 expressions for these quantities in terms of the mode coefficients ca<t, and in Sec. |K 3 
compute the physical energy deposited in a star by an externally applied acceleration a cx t to second order in a ex t- 

Friedman and Schutz |21| derived several general properties of the quantities -E p hy S , >/phys and -E C an,rot- First, the 
physical energy -E p h ys [C] of a perturbation £ ^| has a piece that is linear in £ as well as pieces that are quadratic and 
higher order in 

E phys [C] = E$JC] + E$JC, C] + 0(C 3 ), (Kl) 
and similarly for the physical angular momentum: 

■/physlC] =4hU] + J £U>C]- (K2) 

Second, one has the identity 

E<£ ys [c]^nj^ y M (K3) 

and so the physical energy in the rotating frame, defined as 

-^phys,rot = -^phys ^<^phys C^^v 

has no linear term in £: 

£phy S ,rot [C] = 4hy S ,rot [C, C] + 0(0 (K5) 

4 2 ^ B , rot [c, c] = physic, c] - njMjc a (K6) 

Third, this quantity coincides with the so-called rotating frame canonical energy -E C an,rot[C] defined in Ref. that 
is, the conserved quantity related to the time translation symmetry in the rotating frame by Noether's theorem. 



where 



1. Explicit expressions in terms of £ 

The physical energy of a perturbation can be computed using the procedure of Appendix ^ above, and adding the 
kinetic and potential energies T(t) and V(t). To linear order, the result is 



4hy S [C] = / d 3 xpo(x)(n xx). £ + 2f2 x £ 



= (er,£+B-e) = <£ r ,7r+B.e/2> 



(K7) 
(K8) 



where 



35 This is defined to be the difference between the energy of the perturbed configuration and the energy of the background, 
unperturbed configuration. 
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£ r (x) = fi x x 

is the mode function corresponding to uniform rotation and we have used Eq. (A2). The second order piece is 



(K9) 



E, 



(2) 



1 



1 



1 



physL C,C]- 2\^^/ + 2 ( ^ C - + 2\^ B 'V + 4 (B '^ B ' 
= 5<T,ir) + i(€,C-C) + g(B-€,B.C). 



(K10) 



Here we have used the fact that £ and £ are real to eliminate terms such as ( £ , B • £ ) . The physical angular 
momentum is given by 



Jphys = J dMy x u(y). 



(Kll) 



Here we are using the notation of Append ix where y and u are the inertial-frame location and velocity of a fluid 
element. Using Eqs. ( p7|) and (|j9| ) in Eq. ( Kll ) gives 



U(t) 1 -Jphya= / rf 3 xp (x)(x + |) X fiXX + fJx£ + £ 

and the z-component is given by 

^•^physK] = ^ ' Jphys = O ■ U(t) • Jphys 

= / (fx p (x) [fixx + Ox^]- J7xx + f2x£ + £ 



(K12) 



(K13) 



Expanding this in powers of £ one obtains at linear order using Eq. (K8) the identity (K3). At second order we obtain 



(B • £ , B • 



1 



The physical energy in the rotating frame is therefore, from Eqs. (K6), (K10), and (K14) 

E phys , Iot [C] = \ (t, k) + \ <€ , c • €) + (C 3 ) 

= i (tt , tt) + ~ <£ , c • + i (B • ^ , b • o - i (tt , b • o + o(c 3 ) 



(K14) 

(K15) 
(K16) 



n 



The expression (K15) agrees with the expression Eq. (43) of Rcf. pi[ for rotating- frame canonical energy, as expected 



2. Expressions in terms of mode coefficients 



From the above discussion we can write the total energy of the perturbation, to quadratic order, as the sum of three 
terms, a quadratic rotating- frame energy (physical or canonical), a linear angular momentum term, and a quadratic 
angular momentum term: 



^phys[C] = 4wot[C, C] + fiJ^UC] + ^physIC, C] + 0(C 



(1) 



r(2) 



(K17) 



Of these three terms, the first is diagonalized by the mode basis Caoi the second can be expressed as the coefficient 
CAa of a single Jordan chain mode (corresponding to uniform rotation) , and the third is not diagonalized by the mode 
basis C Aa and contains cross terms between different modes. 



36 Equation (43) of Ref. is actually an expression for inertial-frame canonical energy. The operators which appear in that 
expression, which we will denote by Am, Bin and Cm, are denned such that the linearized equation of motion is Am ■ £j n + 
Bin ■ £ + Cin • £ = 0, where £ in is the inertial-frame Lagrangian displacement. However, one can obtain an expression for 
the rotating-frame canonical energy, as defined by Eq. (64) of Ref. Km, by replacing in Eq. (43) the inertial frame quantities 
Ain, Bin, Cin and £ in by their rotating-frame counterparts used in this paper, A, B, C, and £. Note that in our notational 
convention the operator A is the identity operator. 
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a. Rotating-frame energy 



Consider first the rotating-frame energy. From Eqs. (A5), ( AG), ( A46 ) and (K16) we can write this as 



4wot[CC] = !<C,M.T.C) 



(K18) 



Inserting the mode expansion ( A18 ) and using the definition ( |A14| ) of a right Jordan chain and the definition ( |A5 
of the matrix Ma<j,Bt gives 



E, 



(2) 

phys,rot 



[C, C] = °Aa °Br 



Aa,Br 



-^-Mact,Bt + -K-M.Aa,B{T-l) 



Since the energy is real and M. is Hermitian we can rewrite this as 



E, 



(2) 

phys,rot 



^1 = ! I! C *A* c Bt [(t^A + U B )Ma<t.Bt + iM A a,B{r~l) ~ iM A {<t-X),Bt] 



(K19) 



(K20) 



Ao,Bt 



Since the matrix M. is block diagonal (see Sec. A 5 above), this rotating-frame energy can be written as a sum of three 
contributions, from (i) real-frequency modes, (ii) complex-frequency (unstable) modes, and (ii) non-trivial Jordan 
chain modes. The term corresponding to real frequencies is 



CA\ 



(K21) 



where w e have assumed that one is using a basis which diagonalizes the real-frequencies portion of t he m atrix M_ab 
(cf. Sec. II C 3| above), and ea is the mode rotating-frame-energy at unit amplitude, given from Eqs. (2.36) and ( 2.4(f ) 
by 



E A = UA MaA = UAbA- 

Using the notations ( [2.42 ) we can rewrite the energy ( K2l| ) as 

^ \c a \ 2 uj a b a = 2J \c a \ 2 E a . 



(K22) 



(K23) 



To write down the portion of the expression (K20) corresponding to com plex frequ enci es, we need to introduce some 
more notation. The modes come in pairs (£, uj) and (<?*£*, uj*); see Sees. A 4 and A 5 above. We write A = (aSk), 
where a and 5 label the frequency LUaSk — toas, and k labels the different (degenerate) eigenvectors associated with 
that frequency. The index d can take on the values 1 and —1 such that 



w al — w a 
W a(-1) = ^1 
^>alk ^ak 
£a(-l)fc = (.9*£afc)*- 

It is possible to choose the mode basis such that the matrix M. takes the form 

M a8kM = 6 &b 5 8<-t Skl [ S SA V Hk + Ss,-l?>tk\ ; 



(K24) 
(K25) 
(K26) 
(K27) 



(K28) 



(2) 

cf. Eq. (A5£) above. Inserting this into the expression (K20) gives for the contribution to £p hys rot from complex- 
frequency modes 



L k 



(K29) 



We see that if only a single mode is excited, the corresponding rotating-frame energy vanishes, in agreement with the 
result of Friedman and Schutz pl[ . However, if a mode (£, u>) and also its conjugate pair (g*£*,w*) are both excited, 
then there is a nonzero contribution to the rotating-frame energy. 
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(2) 

Finally, there is a contribution to S phys rot from non-trivial Jordan chain modes. For simplicity, we restrict attention 
here to the case of zero-frequency Jordan chains of length one. For this case the expression ( K2C ) reduces to 



^ Ma\,bo c* A1 c B i, 

AB 



(K30) 



where we have used the identities ( A80 )-( |A82| ). As discussed in Sec. A 5 a above, we can always choose the mode 
basis to diagonalize the matrix 



(3 AB = iM A \,Bo = i (Cai j m ' Cbo) 

= (£ao ! £bo) — (%A1 ' ^ ' 



and the expression ( K30| ) therefore becomes 



(K31) 
(K32) 



(K33) 



where (3 a — Paa are the eigenvalues of the matrix (3. Note that there is no rotating- frame energy associated with the 
cao coefficients. 



To summarize, the total rotating- frame energy is given by the sum of the expressions ( K23 ), ( K29 ), and ( K33 ) 



b. Linear piece of angular momentum 



The linear piece of the z-component of the angular momentum (proportional to the linear piece of the physical, 
inertial-frame energy) can be written as the coefficient of a single Jordan-chain mode. That Jordan-chain mode 
corresponds to the degree of freedom of uniform spin-up Q — ► f2 + AJ7 of the star, and is a zero-frequency Jordan 
chain of length one. As explained in Sec. [D] above, the space of such zero-frequency, length one Jordan-chain modes 



contains all the differential-rotation modes and is infinitely degenerate. In Sec 
compute the basis Xa<t °f left Jordan chains of this space, starting from a specified basis Ca 
Here, however, it turns out to be more convenient to specify the basis by first choosing the basis Xa 
chains, satisfying Eq. (A15), and by defining C Acr to be the dual basis defined by Eq. ( Al6| ). 

The uniform-rotation mode, specified in terms of its left Jordan chain \ a , a = 1,2, is as follows. Let £ r be the 



A 5| above, we explained how to 
of right Jordan chains. 

of left Jordan 



uniform-rotation mode function (|K9|), and define T\ to be any solution to the equation 



c-n =B-e r - 



Define 



and 



Xo = 



-B-£ r /2 



Xi = 



C - B ■ n/2 



(K34) 



(K35) 



(K36) 



Then one can verify using the relation C • £ r = that the left Jordan chain relations • = and • Xi = Xa 
are satisfied. The basis Xao- °f left Jordan chain modes can be completed by adding to the a bove c hain any li nearl y 
independent set of left Jordan chains describing the differential rotation modes. From Eqs. (A19), (K8) and (K35), 
the coefficient c a with a = 1 corresponding to the uniform-rotation mode is 



Cl = (Xo , = «r , 7T + B • |/2) 



E, 



(i) 

phys 



(K37) 
(K38) 
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c. Quadratic piece of angular momentum 



From Eq. (K14), the quadratic piece of the angular momentum can be written as 



(2) 
phys 



[C) C] = (C > K • C) = ^2 C *Aa c BrlCAa,BT, 



A<j,Bt 



where 



K 



-B 
B 



(K39) 



(K40) 



and K-Aa.Br = (Cacf ■• K • C rt)- For non- Jordan chain, real-frequency modes we can write the coefficients /Cab, using 
the orthogonality relation ( 2.32 ) and Eq. ( A32 ), as 



ICab = - (uja+^b) (€a > iB ■ £b) 



i (B.^,B-C B ) 



= ^Ul A b A SAB - \{UA + ^b) 2 (£ A , £i 



- <B • £ A , B • £ B ) 



(K41) 
(K42) 



For regular modes in a rotating star, in the limit O — > 0, the second term in Eq. (K4f ) scales as 0(f2 2 ) and is therefore 
small compared to the fi rst term which is linear in J7. Therefore in general the off-diagonal elements of K-ab are non 
vanishing; see Eq. (3.20) above. 



3. Energy deposited by an external force 



Turn now to the question of how to compute the energy deposited in a rotating star by an externally applied 
acceleration a ext (x, t), to quadratic order in a oxt . The rotating- frame energy can be easily computed, just as one 
computes the energy deposit ed in a non-ro tatin g star. Assuming that there are no unstable mod es, o ne evo lves the 
mode amplitudes using Eq. (2.37) or Eqs. ( A8£ )-( A89), and one computes t he en ergy from Eqs. (K23) and (K33). 

The angular momentum contributions to the total inertial-frame energy (K17) are more difficult to compute. As 
shown above, the quadratic piece of the angular momentum contains cross terms between different modes. Also, in 



evaluating the linear term £ph ys [C] given by Eq. (K38), one needs to solve the equation of motion for the coefficient 
Ci(t) of the uniform-rotation mode to quadratic order in a cxt , i.e., one needs to solve for the back reaction of all the 
other stellar modes on the uniform-rotation mode. The linear equation of motion for the coefficient c\ (t) is given by 
combining Eqs. (A7), (A22), (K35), which gives the equation for angular momentum conservation to linear order: 



ci 



(K43) 



To quadratic order in a cxt (or equivalently in £), one can obtain the equation of motion for c\(t) by using the following 
form of the second-order equation of motion 



£ + B£ + C£ = a< 2 > K, £] + a cxt (x + €,<) + 0(£ 3 ), 



(K44) 



which generalizes Eqs. (2.7) and (|f.l|) above. Replacing the external acceleration on the right hand side of Eq. ( K43| ) 
with the right hand side of Eq. (|K44|) , and using the identity 



( r ,a( 2 )[a]) = -(BiCO 
which is derived below, one obtains the second-order equation of motion for c\ in the form 



ci 



i(B.£,C-0+O(C 3 )- 



(K45) 



(K46) 



The third term on the right hand side here describes the excitation of the uniform-rotation mode by the other stellar 
modes. 
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It is possible to get around these difficulties in the computation of the angular momentum deposited in the star, 
by computing the total torque acting on the star, and by not attempting to compute how the angular momentum 
deposited is distributed between the various stellar modes. The total torque is given byp| 

(K47) 
(K48) 



n J = n j d a xp (x) [x + £(x)] xa oxt (x + £) 

= (£ p , a ext ) + (B • £/2, aert) + (€ P , (£ ■ V)a cxt ) 
Suppose now that the external acceleration can be expanded as 



0(a e a xt,£ 3 )- 



cxt 



e 2 a^ t + 0( £ 3 ), 



and that the response of the stellar modes is similarly expanded as 



(K49) 



(K50) 



The second-order angular momentum deposited can now be computed as follows. First, solve for the linear response 
Cj 1 ] ( t) of the star. Second, compute the angular momentum deposited due to the first term on the right hand side of 
Eq. (K48), which is the standard express ion fo r torque. Third, add to this the angular momentum deposited by the 
remaining two terms in the total torque ( |K48[) , which can be written as 



Aa 



with 



T Aa {t) = I d A x p (x) £ Aff (x) ■ V O • x x a£[(x) 



(K51) 



(K52) 



This procedure yields a expression for the total angular momentum deposited which is a sum over modes (Aa). Note, 
however, that the term Aa in the expression (K51) is not the angular momentum deposited into the mode as 
noted above the total angular momentum contains cross terms between different modes. 
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TABLE I. Coupling coefficient summary for modes of zero-buoyancy stars, to zeroth order in the star's angular velocity. 
Dashes denote coefficients which are forced to vanish by selection rules. Stars denote nonzero coefficients. 





2 pure r modes 


2 regular 


2 hybrid rotational 


pure r mode 






* 


regular 




* 


* 


hybrid rotational 




* 


* 
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